text:
# 1 Introduction

## 1 Introduction

Lattice-QCD theory predicts the disappearance of the hadronic phase of matter once the energy density exceeds a critical value of about 1 GeV per fm3 , giving rise to a continuous, deconfined QCD state that is governed by the elementary interaction of quarks and gluons. To recreate this phase in the laboratory one collides heavy nuclei at relativistic energy with the goal of ascertaining the QCD predictions, and to pin down the decay point from the partonic to the hadronic phase by obtaining estimates for the transition temperature and energy density. The CERN SPS Lead (208Pb) beam facility provides for a top energy of 158 GeV per projectile nucleon, corresponding to a total internal CM-system energy of about 3.5 TeV, to heat and compress the primordial reaction volume. A tenfold higher $`\sqrt{s}`$ is reached at the BNL heavy ion collider RHIC. In fact calorimetric data  show that the average transverse energy density exceeds about 2.5 GeV/fm3 already at the SPS in central Pb+Pb collisions. Moreover the study of $`J/\mathrm{\Psi }`$ production  demonstrates a suppression of the yield in such collisions, characteristic of the QCD ”Debye” screening mechanism expected in a deconfined partonic medium .

This lecture will deal with physics observables that could tell us about that medium. As an introduction I will look at the reaction dynamics in AA central collisions and then turn to my main topic, bulk hadron production and the location of the QCD phase boundary. We shall see that the production rates of all strange hadronic species play a prominent role, as a signal that helps to understand dense QCD matter, and its decay.

### 1.1 Fireball Matter Dynamics

The initial Woods-Saxon nucleon density distributions, with average energy density of about 0.16 GeV per cubic fermi, impinge onto each other leading, at first, to nucleon-nucleon collisions occuring concurrently at the microscopic level as enveloped by the overall impact geometry of the target/projectile nuclei. This first generation of binary nucleon-nucleon encounters will involve essentially all the incident 2 A nucleons if we consider head-on collisions of nuclei with mass number A. In Pb+Pb collisions we thus have about 200 primary such nucleon-nucleon collisions of first generation. If we could now somehow stop the reaction dynamics, letting the reaction products escape to a detector system that identifies them, we would expect to register a final multiparticle state with a composition closely resembling A-times the (well known) outcome of a nucleon-nucleon collision at similar center-of-mass energy $`\sqrt{s}`$: nothing new!

However what makes relativistic nuclear collisions interesting is the fact that, in reality, the first generation set of A microscopic nucleon-nucleon interaction systems will immediately re-interact while still in statu nascendi of its ”pending” output of asymptotically distinguishable, ”on shell” reaction hadrons which would consist of about 10-20 produced mesons and baryon-antibaryon pairs. However, before being fully formed they run into a secondary generation cycle of subsequent collisions within the nuclear density distributions of the heavy nuclei. In fact there may be up to 6 secondary, subsequent collision generations of such pre-exited (not finally formed) microscopic collision volumes in a central Pb+Pb collision. Moreover they will occur in successive time steps spaced by fractions of a fm/c only \- due to the relativistic Lorenz contraction of the nuclear density profiles of the Pb nuclei in beam direction ($`\gamma _{CM}=9`$ at the SPS).

Something new! As we have no experimental knowledge of such secondary collisions of partly incompletely formed (”off shell”), partly decomposed hadrons we have to withstand the temptation to capture the overall reaction dynamics in a classical billiard ball cascade of subsequent generations of an inelastically multiplying gas of known hadrons. Into such a model we would insert the known, vacuum elementary cross sections at each microscopic binary encounter and proceed via Monte Carlo probability sampling methods. Such approaches are called microscopic hadron transport models . From the above we would be surprised if they could give a correct description of the final outcome of a highly relativistic heavy ion collision because their set of microscopic degrees of freedom (isolated binary collisions of ”on shell” hadrons) fails to capture the unknown nature of secondary, tertiary etc. encounters in the dense medium. These encounters, moreover, occur unresolved in time which should lead to quantum mechanical coherence, which might render the entire picture of isolated, sequential microscopic cascade-sub-processes obsolete. Let us, therefore, conclude that the overall large interaction volume of a central Pb+Pb collision will (after about 2-3 fm/c of interpenetration and reaction time have elapsed) be composed of a hitherto unknown state of strongly interacting matter that, however, contains all the quantum numbers and relative center of mass energy of the initial nuclear projectiles. We thus suspect that it will, quite generally, feature a high spatial density and, similarly, a high energy density. We therefore sometimes call this short-lived object a ”fireball”. The state of strongly interacting matter inside it is the topic of modern QCD theory  and of our experiments.

Now there is something simple about this picture of dense matter in a fireball that experiments can check quickly. A new ”state” of ”matter” supposedly has been formed by fusing spheres of target and projectile cold nuclear matter which are initially located at opposite ends of longitudinal momentum space in the CM-system. By the symmetry of Pb+Pb collisions, a fused fireball should occupy a common momentum space volume centered in the CM-system at zero longitudinal and transverse momentum. Fig.1 shows the distribution in longitudinal phase space (measured here by the rapidity variable $`y=0.5ln[(1+\beta _L)/(1\beta _L)]`$ for the negative and neutral K-mesons, and for the anti-strange hyperon $`\overline{\mathrm{\Lambda }}`$ : all are Gaussians well centered at $`y_{CM}=0`$. Note that the $`K^{}`$ and $`\overline{\mathrm{\Lambda }}`$ consist only of newly created quarks that were not brought into the fireball as initial nucleon valence quarks. For an ideal, spherically symmetric fireball we would also require an isotropic momentum distribution. Looking simultaneously at the transverse momentum (not shown here) and at the y-distributions of these particles we see that this is not strictly the case (the fireball is longitudinally stretched and looks like a fire-football of US vintage) - but closely enough .

Furthermore we can experimentally gather almost all hadrons created in a central Pb+Pb collision at top SPS energy: it is about 2500 of them! Superimposing the total transverse energy in the football, carried by the various hadronic species we can estimate  its total average transverse energy density in space. This is the 2.5 GeV/fm3 that we mentioned initially . The fireball should thus initially exceed the QCD phase boundary between hadrons and quark-gluon matter, which is in the vicinity of 1 GeV/fm3! Then it will expand explosively, back to hadrons via the phase transformation.

### 1.2 Deconfinement and Phase Boundary: Signals

How do we get experimental signals that elucidate the state of matter in the fireball maximum density stage? Like with supernova analysis, we may try two ways: observe the primordial ”light curve” (i.e. electromagnetic and neutrino radiation) or look at the bulk explosion material (e.g. expansion modes, element composition etc.). For the nuclear fireball the neutrino signal is out but we are left with ”black body radiation”: directly emitted photons and lepton pairs which leave the fireball while it is hot and dense. The corresponding experiments  are very demanding but successful. However I will not cover them here for lack of space.

Focusing on experiments that analyze the fireball material after expansion and cooling I briefly mention the suppression of the J/$`\mathrm{\Psi }`$ yield confirmed experimentally  to occur in central Pb+Pb collisions at top SPS energy, 158 GeV per nucleon. Still discussed vigorously in the community, we may take such observations to be, at least, consistent with creation of a fireball energy density well in excess of 1 GeV/fm3 at which QCD predicts deconfinement thus dissolving the $`J/\mathrm{\Psi }`$. This observation agrees with our above estimate of the transverse energy density as derived from the total final phase space density of the produced hadrons.

If one tentatively takes for granted such indications of a deconfinement state at top SPS energy one expects, likewise, to receive signals of the parton to hadron phase transformation bound to occur once the primordial high density state expands and cools toward the critical temperature and energy density. This brings me to the principal topic of the following chapters, in which I shall develop this observable in broader detail: could it be imagined that the detailed composition of the fireball decay products, in terms of abundances of the various hadronic species, captures the parton to hadron QCD phase transformation period of the dynamical fireball evolution? Or, more generally speaking: we expect to create, at least, a novel state of strongly interacting matter in the short lived fireball which must then decay. We infer from the above considerations of microscopic fireball dynamics that this high density state should feature novel degrees of freedom at the microscopic level. These may be partons (at top SPS energy and beyond), or shadows of the familiar hadron spectrum as modified by interaction with their dense, surrounding fireball medium (as might occur at lower collisions energies), and by the microscopic reaction cycles during the interpenetration stage of the nuclear density distributions. In any case it is the decay of this unknown coherent quantum mechanical state to a quasi classical state with familiar degrees of freedom, in our case a hadron and resonance gas, that we have to consider. The composition of that hadron gas is detectable: our signal.

Let me make a few remarks referring to QCD folklore. Almost trivially, all conceivable states of strongly interacting matter are falling under QCD governance. Thus in a decay of a state higher in energy density and temperature, upon fireball expansion, to the lower density hadron gas it is just the set of proper QCD degrees of freedom that is changing - we are thus dealing with a QCD phase transformation. At the hadron gas level, QCD resides in the spectrum of mass, spin, isospin, flavour, width etc. of the vast array of hadronic species. Thus in a QCD parton to hadron phase transformation, an initial fireball ensemble of flavoured plus coloured quarks and antiquarks, and of colour-anticolour carrying gluons is, sloppily speaking, ”looking down” at the QCD realization below: the hadron spectrum. The ensuing, colour neutralizing ”condensation” process of partons leads to hadronic objects featuring a spectrum of specific flavour-colour-spin-momentum internal compositions that absorb the initial partonic degrees of freedom. Physics experience suggests that such a process, occuring in an extended volume, should be governed by statistics. In the partonic view ”from above” a light ($`u\overline{d}`$ positive) pion will be certainly the easiest way to hadronize in a predominant $`u,\overline{u},d,\overline{d}`$ population (likewise for the other pions). In strong contrast, a heavy ($`sss`$) Omega hyperon will be a highly unlikely hadronization outcome, in view of competing $`K^{}`$ and $`K^0`$ mesons that could take care of the s-quarks combining them with the more abundant $`\overline{u}`$ and $`\overline{d}`$ quarks. Net result of all of this: the hadronization transition populates the hadronic spectrum in order of the relative statistical weight of the hadrons. We call such a process ”phase space (statistical weight) dominated”. It will thus create a hadron gas ensemble of maximum entropy, a decoherent, classical system  which exhibits a characteristic ordering pattern concerning the relative abundance of each hadronic species. It is, thus, not surprising that the abundance spectrum (expressed in terms of the various average hadronic multiplicities per collision event) obeys a statistical Gibbs ensemble , either in the so-called ”canonical” or in the ”grand-canonical” realization.

The statistical analysis of the composition of multihadronic final states was pioneered by Hagedorn in the 1960’s . It was revived in the last decade, applying it to the hadron populations observed in relativistic nucleus-nucleus collisions as well as in $`pp`$ and $`e^+e^{}`$ annihilation . In this approach one captures the temperature and energy density prevailing at birth of the multihadronic final state, i.e. the point in the fireball dynamics where it decouples, by decoherence, from the novel state of high energy density/temperature created in the early phase of the dynamics. This is clearly an interesting signal!

## 2 Hadron Multiplicity and Strangeness Enhancement

The first SPS experiments with 32S-beams at 200 GeV/A showed an enhancement of various strange particle multiplicities, chiefly $`K^+,\mathrm{\Lambda }`$ and $`\overline{\mathrm{\Lambda }}`$, relative to pion multiplicities, in going from peripheral to central S + (S, Ag, Au) collisions . This observation appeared to be in-line with the pioneering analysis of Rafelski and Müller  who first linked strangeness enhancement to the advent of transition from the hadronic to a partonic phase. This offered lower effective $`s\overline{s}`$ threshold, shorter dynamical relaxation time toward flavour equilibrium, and an additional, nontrivial effect of relatively high net baryon number or baryochemical potential: the light quark Fermi energy levels move up, perhaps even to the $`s`$-quark mass at high $`\mu _B`$, and the Boltzmann penalty factor for the higher mass $`s\overline{s}`$ pair creation might be removed. This latter aspect was mostly ignored in the late 1980s but receives fresh significance as we become increasingly aware of the crucial role of $`\mu _B`$.

In this section, and in sections 3 and 4, I will present a sketch of our recent progress, both in gathering far superior data and in the understanding of the statistical model that was rudimentarily anticipated in such early strangeness enhancement speculations.

Fig.2 shows a modern version of the typical strangeness enhancement phenomena. $`\varphi `$-meson to positive pion multiplicity ratios obtained by NA49  in Pb+Pb SPS collisions at 158 GeV/A (corresponding to $`\sqrt{s}`$=17.3 GeV) are plotted for a sequence of collision centrality conditions from peripheral to central. At the peripheral end the minimum bias p+p point  matches with the trend. The centrality scale employed here is, at first, the number of participating nucleons (left hand side of Fig.2). The raw data centrality bins are ordered in NA49 data by decreasing projectile spectator energy as recorded in a zero degree calorimeter. This information is converted to mean participant nucleon number by a VENUS calculation from which one also derives the mean impact parameter b, corresponding to each centrality bin. Neither of these scales turn out to be satisfactory  in merging data from central light nuclei collisions such as C+C with the various centralities of the Pb+Pb collisions. For example, a central C+C collision has $`b2`$ fm and $`N_{part}18`$ but on a b scale the $`\varphi /\pi ^+`$ value is about 50% lower than the $`b=2`$ result for the much heavier Pb+Pb system. Inversely on the $`N_{part}`$ scale: $`N_{part}`$=18 corresponds to very peripheral Pb+Pb and the central C+C result is about 50% higher than the Pb+Pb curve. A central collision of a relatively light nuclear pair thus behaves quite differently from a very peripheral heavy nuclear collision where only the dilute Woods-Saxon density tails interact! The scale of the right hand side of Fig.2 is an intuitive guess  to represent the relative compactness, or volume-to-surface ratio of the primordial interaction volume, by the variable R-$`b`$/2, where R is the radius of the colliding nuclear species. It might be connected with the energy or collision density reached in the primordial collision volume (see chapter 5). We see that the central light nuclear collision data now merge with the Pb+Pb centrality scale.

The ”strangeness enhancement factor” is given oftentimes as the production ratio of AA central/(pp min. bias times 0.5 $`N_{part}`$). In the case of Fig.2 it would be roughly 2.5. We also have a systematic study of strangeness enhancement in Au+Au collisions at the BNL AGS energy of 11 GeV per nucleon where this factor is about three . Multistrange hyperons  show factors between 4 and 15 at top SPS energy.

Bulk strangeness enhancement in central collisions is a nuclear feature, absent in pp collisions. Of course we lack a detailed picture about ”centrality” in pp collisions but we could still employ e.g. the total charged particle multiplicity to select more or less ”violent” collisions. Fig.3 shows the $`K^+/\pi `$ ratio of pp at 158 GeV versus charged particle multiplicity to be essentially flat . Similar findings are made up to Tevatron $`\sqrt{s}=1.8`$ TeV $`p\overline{p}`$ collisions: the $`K^+/\pi `$ ratio is 50% higher here but also almost independent of $`N_{ch}`$ .

A picture emerges in which strangeness enhancement, or more generally speaking the yield order in the overall bulk hadron population is connected with ”sequentiality” of interactions at the microscopic level, i.e. with the fact that a number of successive collisions occurs in narrow intervals of time if one may employ a naive Glauber picture. This number reflects the size and density of the primordial interaction zone. It thus seems that a certain (extended) size and a certain interaction rate (and thus energy density), typical of central nuclear collisions, conspire in strangeness enhancement. Unfortunately we do not really know how to describe a second, third etc. collision of a hadron, within fractions of a fm/c space-time distance. We would then have probably resolved the key issue: does it dissolve into a parton cascade from which the final hadrons are reconstituted? Proton-nucleus collisions should hold a key to this question but nobody has succeeded in isolating the second, third, n’th successive collision of the projectile, as of yet .

At the moment we thus forgo pA as an intermediate step although it certainly also features changes in the hadronic production ratios  and base the analysis on comparing pp to AA. Fig.4 shows the hadronic multiplicities, from pion to cascade hyperon, obtained by NA49 for min. bias pp at $`\sqrt{s}`$=17.3 GeV . The data are confronted with the Hagedorn statistical model in its canonical Gibbs ensemble form as employed by F. Becattini , leading to very good agreement (as it is well known also for other elementary collisions and energies ). The three parameters are T=186$`\pm `$7 MeV, a reaction volume of about 7 fm3, and a total of about 0.5 $`s\overline{s}`$ pairs. The apparent validity of a statistical weight-dominated picture of phase-space filling has been considered a puzzle already since Hagedorn’s time. It is clear, however, that the apparent canonical ”hadrochemical equilibrium” pattern can not result from ”rescattering” of produced hadrons: there is none. In Hagedorn’s view  a creation ”from above” must hold the key to the apparent maximum entropy state, i.e. the QCD process of hadronization . In it a medium of unknown composition (the ”fireball” or ”string” or ”cluster”) undergoes a quantum mechanical decay which is dominated by the phase space weight distribution that corresponds to the spectrum of hadrons and resonances . The canonical abundance pattern and the T-value are a fingerprint of QCD hadronization - do AA collision data at high $`\sqrt{s}`$ also confirm this picture (they must, of course, also result from a hadronization process)?

## 3 AA collisions in the Grand Canonical Model

Fig.5 shows the grand canonical fit by Becattini to the NA49 data from central Pb+Pb at 158 GeV/A . The temperature is 160 $`\pm `$ 5 MeV and $`\mu _B`$ =240 MeV; besides, this model employs the much discussed strangeness undersaturation factor $`\gamma _s`$=0.8.

Leaving the second order concern about $`\gamma _s`$ to the theoretical community I note here that Braun-Munzinger et al.  fit a set of data at the same top SPS energy without introducing a $`\gamma _s`$; they report T=170 $`\pm `$5 MeV, at $`\mu _B`$=270 MeV, close enough. There are also studies of the new RHIC STAR data  at $`\sqrt{s}`$=130 GeV by this model  and by Kaneta and Xu , averaging at 175 $`\pm `$ 5 MeV and $`\mu _B`$=48 MeV. And the new data of NA49  at 80 and 40 GeV/A have resulted in fit values of T=155 MeV, $`\mu _B`$=270 MeV and T=150 MeV, $`\mu _B`$=395 MeV, respectively . I will return shortly to a further discussion of the grand canonical approach but wish to, first of all, show an overall impression from these analyses which are confronted in Fig.6 with the new lattice QCD calculations at finite $`\mu _B`$ by Fodor and Katz .

The latter predict the T, $`\mu `$ dependence of the QCD phase transformation which in this model consists of a crossover for all $`\mu <`$ 650 MeV, i.e. in the SPS-RHIC domain. Note that physics observables can change rapidly in a crossover, too: the familiar steep rise of e.g. lattice $`ϵ/T^4`$ at $`T_c`$ does not, by itself, reveal the order of the phase transformation . Anyhow: the hadronization points from grand canonical ensemble analysis merge with the phase transformation site of lattice QCD at top SPS and RHIC energy. Quite a remarkable result, but also a plausible one  if we recall that Ellis and Geiger did already point out in 1996 that hadronic phase space weight dominance appears to result from the colour-flavour-spin-momentum ”coalescence” of partons that occurs at hadronization . Unfortunately a rigorous QCD treatment of the parton to hadron transition is still missing.

At this point the following objection is always raised: if the same basical model describes hadronic yield ratios in $`pp`$, $`e^+e^{}`$ and in central AA collisions, Figs. 4 and 5, what is special about AA, as you will not tell us now that a QGP is also formed in pp?! Answer: on the one hand both collision systems reveal the QCD hadronization process which features, furthermore, the Hagedorn limiting hadronic temperature $`T_H`$. At top SPS and at RHIC energy $`T`$ (hadronic ensemble) $`T_HT_c`$ (QCD), this is the common feature; it should not be a coincidence. On the other hand hadronization appears to occur under distincly different conditions in AA collisions, as captured in the transition from a canonical to a grand canonical description. Inspection of Fig.4 and 5 shows that the hadronic population ratios are quite different: the falloff from pions to strangeness-two cascade hyperons in the former case is about four orders of magnitude whereas it reduces to three in the grand canonical situation: strangeness enhancement! In the canonical case the small reaction fireball volume is strongly constrained by local conservation of baryon number, strangeness neutrality and isospin whereas these constraints fade away in the grand canonical ensemble which represents a situation in which, remarkably, these conservation laws act only on the average, over a rather large volume, as captured by a collective chemical potential $`\mu _B`$. This leaves one global quantity $`\mu _B`$ essentially in charge of all the conservation tasks. Note that the statistical model does predict nothing, it merely captures this most remarkable feature of the hadron gas emerging after hadronic freeze-out. Its observed success implies some kind of long range collective behaviour in the hadronizing source, the origin of which is yet unknown, but must be specific to central AA collisions. We shall address this topic in the following Sections but note, for now, that Strangeness enhancement is the fading away of small volume canonical constraints, in the terminology of the statistical model .

This aspect has been recently studied in all detail by Cleymans, Redlich, Tounsi and collaborators . Fig.7 illustrates their results concerning the transition from canonical to grand canonical behaviour with increasing number of participants, i.e. overall ”source” size. It is intuitively clear that it should occur, first, in singly strange hadrons, the increase occuring with offset (but having a larger specific effect on the yields per participant) in S=2,3 hyperons.

A further, appropriate critical question: how can we understand the other aspect of Fig.6, i.e. the steep falloff from the QCD transition domain occuring at the lower SPS energies? We even have a further GC analysis, at top AGS energy, by Stachel , for central Si+Au collisions at 14.6 GeV/A, shown in Fig.8.

The result is $`T`$=125 MeV, $`\mu _B`$=540 MeV, far below the $`T`$ scale of Fig.6. The picture of a direct parton to hadron transition is intuitively inapplicable at these lower energies. Still the overall dynamical trajectory that ends in hadronic chemical (abundance) freezeout should arrive there ”from above” as hydrodynamical models  show. The hadronic system is initially very dense along its trajectory, it is not a hadron gas, and it is thus a quantum mechanical coherent state composed of excited and in-medium modified hadrons that decays to the finally observed classical hadron gas ensemble. This decay might also ignore the classical concept of a relaxation time. Recall the nucleus, also still a dense system: we do not invoke relaxation time in a transition within such a quantal medium, such as $`\beta `$-decay. And yet ”Fermis Golden Rule” asserts that the transition strength depends ”only” on the squared matrix element times final state phase space volume weight plus global conservation laws. And we know that the phase space factor oftentimes far overrides the matrix element, in the net decay strength. High density hadronic matter behaviour is essentially unknown: an old and new research paradigma but at present we state that it is not a classical hadron-resonance gas and that it may decay quantum mechanically, the decay products well described by a grand canonical ensemble of hadrons and resonances. At top SPS and at RHIC energy, in turn, the increasingly ”explosive” nature of partonic and hadronic expansion may almost instantaneously dilute the hadronizing source toward chemical freezeout, as indicated by T(GC)$``$ T (Hagedorn) $`T_c`$ (QCD). We may, thus, directly look here at the QCD parton to hadron phase boundary, located at T $``$ 170 MeV.

## 4 Working with the Grand Canonical Ensemble

### 4.1 The Hagedorn Fireball

The previous sections have described our recent attempts to understand the order of multiplicity per collision in which the various hadronic species are populated. We have referred to the statistical Hagedorn model  in its modern versions. It has turned out that elementary collisions, such as $`pp`$, $`p\overline{p}`$, $`e^+e^{}`$ $``$ hadrons, are well described in the canonical version of the statistical model, whereas central nucleus-nucleus collision hadronic final states appear to obey grand canonical statistical order. Let me emphasize again that in this model we consider the finally observed multihadron state to result from the decay of a quantum mechanically coherent fireball stage that resides ”above” the finally observed, ”frozen-out” classical hadron gas, in terms of energy density and temperature. Its decay occurs at a certain, late stage in the overall dynamical evolution, imbedded into an overall cycle of initial interpenetration, fireball stage of maximum energy density, and expansion dilution and cooling. The composition and structure of matter in the intermediate fireball stage is the object of foremost interest. Its direct radiation output signals form one set of appropriate observable signals . The other signals are derived from the fireball freeze-out decay into hadrons.

The freeze-out product of the fireball is sufficiently dilute to be quantum mechanically incoherent, thus being describable as a classical hadron gas. The multiplicities, and multiplicity ratios of the various hadronic species in this gas form an image of the instant of decay: we thus study the decay properties of the unknown state of matter in the high density fireball: the conditions prevailing ”at birth” of the frozen-out hadron gas, common to all hadrons emerging from the fireball. These common conditions are, chiefly, the temperature, energy density and net baryon number density. They are captured by the statistic ensemble analysis. Its success in describing the composition of the frozen out hadron gas (residing in the multiplicity ratios of dozens of hadronic species, from pions to multiply strange heavy hyperons) shows that the entire, at first sight exceedingly complicated final state, of up to several thousands of hadrons, can be understood from a few common macroscopic parameters. To summarize: it is the underlying assumption of Hagedorn analysis that the ”decay from above” is dominated primarily by the energy density (and conserved quantum numbers) at the instant of decay, and by the statistical weight distribution offered by the hadronic spectrum into which the decay occurs. Immediately after decay the hadronic gas thus exhibits a statistical equilibrium population, as captured by a Gibbs ensemble of mixed hadronic species. At the risk of overdue repetition, I stress again that this apparent ”thermal” equilibrium is a result of the decay process, the nature of which lies well beyond the statistical model which ”merely” captures the apparent statistical order, prevailing right after decay. The hadrons and their weight distribution are ”born into” this state ! The observed equilibrium is, thus, not achieved by inelastic transmutation of the various hadronic species densities, in final hadron gas rescattering cascades, i.e. not by hadron rescattering approaching a dynamical equilibrium. The system explodes and cools rapidly after the initial hadron gas formation phase , the canonical or grand canonical order staying frozen-in throughout subsequent expansion (while the momentum spectra etc. still get modified by elastic interactions, resonance decays etc.). Hadro-chemical composition freeze-out thus occurs prior to the final decoupling from all strong interaction (spectral freeze-out). In a typical central Pb+Pb collision at top SPS energy, 158 GeV per projectile nucleon, hadro-chemical freeze-out occurs at a temperature of $`165\pm 5`$ MeV (thus capturing the hypothetical QCD phase boundary) but final, spectral decoupling occurs at about 110 MeV , as we learn from hadronic spectra and correlations.

I have inserted this somewhat lengthy section to avoid the misunderstandings and controversies that have accompanied the Hagedorn model ever since its inception. In the 1970’s particle physicists were seeking an understanding of proton collisions in analogy to the elementary Feynman graphs of QED, thus trying to insolate similarly ”elementary” hadronic processes, connecting the initial and final constituent (massive) quarks. However the most probable outcome of a $`pp`$ collision, namely to go to a ten-hadron final state (of statistically varying microscopic composition) clearly defied such a picture, thus being called ”background” reactions, outside of the primary research focus. On the other hand Fermi, Landau and Hagedorn just put the emphasis on this side of hadron collisions, guided by the intuition that an overall process of ever increasing density of potentially coexisting microscopic subprocesses should approach the statistical situation: all subprocesses would feed into a symbolic, intermediate ”compound state” (reminiscent of Bohrs excited nuclear compound state)that was called fireball. Its decay would then feed into a frozen-out, statistically ordered hadron gas. The appropriate formulation turned out to be a canonical Gibbs ensemble . We have seen the modern version of this model in Fig.4.

Furthermore, the Hagedorn concept of freeze-out from the symbolic fireball (of microscopically unknown degrees of freedom) to a classical hadron gas of statistical weight dominated composition - that then travels on to final observation unchanged by the subsequent dynamical evolution - took a long time to be comprehended. On the other hand we are cordially familiar with the same picture from explosive Big Bang nucleosynthesis. We know that the cosmic average proton to helium composition ratio froze out to its observed ratio once the inelastic transmutations in the cosmic fireball (among the various light nuclear species) stopped at about 1 MeV temperature - thus presenting to us highly relevant data concerning this dynamical stage. With further expansion this chemical composition travels on unchanged (frozen-in), while the spectral temperature of the cosmic inventory has dropped down to a few $`10^4`$ eV. Likewise, hadro-chemical fireball freeze-out creates hadron composition ratios that travel onward throughout further hadronic expansion.

### 4.2 The Grand Canonical Ensemble

Let us now turn to the formalities of grand canonical ensemble analysis. One starts from the formulation of a partition function which specifies the relative weight $`Z_i`$ (it counts the sum of possible states) for each particle (or resonance) species $`i`$ in a multihadronic mixed gas at temperature T:

$$lnZ_i=\frac{g_iV}{6\pi ^2T}_0^{\mathrm{}}\frac{k^4dk}{E_i(k)[exp\{(E_i(k)\mu _i)/T\}\pm 1]}$$
(1)
where $`g_i`$ is the statistical Land$`\stackrel{´}{\mathrm{e}}`$ factor (the ”statistical degeneracy”) of species i, V is the total common volume shared by all species, $`E_i^2(k)=k^2+m_i^2`$ the total energy square of species i at momentum k, and $`\mu _i=\mu _BB_i+\mu _sS_i+\mu _II_i`$ the ”chemical potential” of species i. The latter is the typical, unfamiliar ingredient of the grand canonical ensemble. Please consult textbooks describing its occurence in $`Z_i`$ as a result of adding a so called Lagrange multiplier to the Lagrange density of the system in order to enforce global conservation of certain net quantum numbers that are specific to the total system contained in the total volume V. For our case of ”hadrochemistry”, these are the net quantum numbers initially carried into the fireball volume by the incident nuclei. Their total baryon number B, total strangeness S, and total isospin (Z-component) I are initially given by the total participating nucleon (=baryon) number, S being zero (nuclei have zero strangeness), and by the net isospin of all participant nucleons: I=(Z-N)/2. These net quantum numbers are conserved by strong interaction, and accompany the collision volume throughout the fireball evolution, during which they will be re-distributed over the entire ensemble of the hadronic species that is being created. Thus the three components of $`\mu _i`$ in Zi represent the net impact of the overall quantum number conservation on each separate species i. We notice, most significantly, that the conservation laws are thus enforced on average only, not locally but over the entire fireball volume. This total volume thus acts somehow coherently, balancing globally the quantum number exchanges that have occured microscopically during fireball dynamics preceding the decoupled hadron gas that is described by this set of Zi’s. From our considerations above, and following Hagedorn’s advice , we state again that this global coherence is not achieved by inelastic relaxation toward equilibrium within the frozen-out hadron gas ensemble but must be a characteristic feature either of the quantum mechanically coherent high density fireball preceding freeze-out, or by the mechanism of its decay process to the decoupled hadron gas (or by both acting together). After decay the resulting hadron gas is imprinted with these global coherence aspects of its birth process. The grand canonical ensemble description merely ”takes note” of the freeze-out product of such a process, analyzing a temporary snapshot of the apparent order in the frozen-out state, ignorant of its dynamical origin.

A short comment concerning the terminology: chemical ”potential”. You see from equ. (1) that the exponential ”penalty factor” in the denominator, that expresses the cost of realizing a certain hadron with total energy Ei(k) within a bath of temperature T (characterizing the average kinetic energy, or energy density), gets modified by $`\mu _i`$. This modification thus takes into account the ”affinity” - to employ the classical terminology - that the medium offers to species i. In modern terms we realize that $`\mu _i`$ in the exponent acts like an average potential for species i in the medium, modifying its vacuum energy $`E_i(k)`$, to now read $`E_i(k)\mu _i`$. The medium thus modifies the cost of realizing species $`i`$ at momentum $`k`$. An ”affine” medium (positive $`\mu _i`$) reduces the penalty factor and thus increases the relative weight $`Z_i`$.

From the partition function in (1) we derive the distribution of number densities of the species i in the medium,

$$n_i=\frac{g_i}{2\pi ^2}\frac{k^2dk}{exp\{(E_i(k)\mu _i)/T\}\pm 1}$$
(2)
The contribution of species i to the overall energy density $`ϵ_i`$ is obtained by inserting the factor $`E_i(k)`$ in the numerator of the integral. From the number density of species i we get its total multiplicity per event by multiplying $`n_i`$ with the total reaction volume V.

Now we want to analyze experimental data, with the set of $`N_i`$ hadronic multiplicities supposed to be known. From this we have to determine the basic parameters V, T and $`\mu _B`$, $`\mu _S`$ and $`\mu _I`$. It turns out that we can, before proceeding to a fit to equ. (2), reduce the three chemical potentials to one remaining quantity. We do this by exploiting global baryon, strangeness and isospin conservation in a set of coupled equations which, in the end, express $`\mu _s`$ and $`\mu _I`$ in terms of $`\mu _B`$ which thus remains as the only independent potential, expressing all conservation tasks. Three parameters remain: V, T and $`\mu _B`$. We have shown in Fig.5 the outcome of such a fit procedure as applied to a set of NA49 hadron multiplicity data for central Pb+Pb collisions at top SPS energy. Repeating at various energies  we get the systematics of hadronic chemical freeze-out conditions as represented in the T, $`\mu _B`$ plane (Fig.6). That is, basically, how it works.

An easy qualitative understanding of the procedure that determines the two parameters T and $`\mu _B`$, in a first cycle (which then leads to fixing V in a second step) is obtained from considering the special case of particle to antiparticle multiplicity ratios. For easy approximation we go from Bose or Fermi statistics (resulting in the $`\pm 1`$ term in the denominator of equ. (2)) to Boltzmann classical statistics dropping that term. We consider the ratio of a certain particle-antiparticle multiplicity,

$$<n_i>/<\overline{n}_i>=n_iV/\overline{n}_iV=n_i/\overline{n}_i$$
(3)
in which the common total volume parameter drops out. Consideration of the appropriate ratio of two integrals (2) shows that almost everything drops out because statistics and phase space are identical for particle-antiparticle pairs. What remains is

$$<n_i>/<\overline{n}_i>=exp\{(\mu _i\overline{\mu }_i)/T\}=exp(2\mu _i/T).$$
(4)
If we do this for several such ratios in combination, say $`K^+/K^{},p/\overline{p},\mathrm{\Lambda }/\overline{\mathrm{\Lambda }},\mathrm{\Omega }/\overline{\mathrm{\Omega }}`$ etc., we get a set of equations from which to determine the maximum likelyhood $`\{\mu _B,T\}`$ combination. In a second step we fix the volume parameter V by fitting equ. (2) to the pion multiplicity $`<n_\pi >`$ with known $`\mu _B`$ and $`T`$.

In the real procedure one has to go through one additional complication. The finally observed hadrons are not really the ones that existed at freeze-out, T=165 MeV or so, in the created hadronic gas. It is composed, at first, of excited hadronic states such as $`N^{},\mathrm{\Delta },K^{}`$ etc. and of resonances like $`\varrho `$ and $`\omega `$. It is their population that obeys the set of equations (2)! After being established at freeze-out this population decays in well known ways to the finally observed $`\pi ,K,\mathrm{\Lambda },p,\overline{p}`$ etc.. These latter particles have never been in a state that obeys the grand canonical multiplicity order! They serve as an observational input to a grand canonical ensemble fit via an attached procedure that was invented by Wroblewski , which relates the final, observed set of hadron multiplicities to a set of excited hadron and resonance multiplicities that forms the initial hadron composition at freeze-out, via the known decay branching ratios. This important backward-transformation was oftentimes ignored in early applications of the statistical model, causing considerable confusion.

## 5 The Onset of Grand Canonical Order in AA Collisions

We infer from the previous section that grand canonical ensemble analysis captures a snapshot at birth of the classical hadron gas, emerging from central collisions of heavy nuclei. The apparent most interesting feature was the quantum number conservation on average only, over a large volume of microscopic ingredients. We have claimed that this coherence is not a property of the frozen-out hadron gas but should reside in its birth process ”from above”, referring hypothetically, to the dynamics of fireball evolution preceding decoupling and freeze-out. We have, thus, invoked a prior dynamics: coherence in a large volume, high energy density fireball (of unknown matter composition) which, upon decay, creates grand canonical order. The grand canonical analysis merely takes note of this order; it does not know about our additional, hypothetical speculation concerning its dynamical origin.

On the other hand, as we are interested, in particular, in just this dynamical origin (as it might tell us about the properties of the unknown fireball matter above: its decay temperature and chemical potential), we would now want to inspect a wider variety of fireball hadronization processes. We thus turn to elementary hadron collisions, and to $`e^+e^{}`$-annihilation, in which the interaction volume is very small in comparison to AA fireball volumes. Nevertheless, hadronization must occur under similar QCD dominance here. This leads us to expect that the statistical temperature, which is the decay image of the critical energy density at which any QCD hadronization occurs, must be common to both elementary and AA collisions at $`\sqrt{s}20`$ GeV (top SPS to RHIC energy). Because we assume, from all existing observables, that we are witnessing a parton to hadron QCD phase transformation at these high energies. Indeed the temperature parameters all fall within $`170\pm 10`$ MeV, from $`e^+e^{}`$ to Au+ Au .

The apparent difference between small and large freeze-out volumes appears, thus, to arise not from the decay temperature nor from the final hadronic level spectrum but from a specific difference in the size of the coherently decaying system. Indeed, in the canonical decay mode the ”volume” amounts to a few fm3 only. The relevant conserved quantities are forced in the canonical case to be conserved locally - if an Omega hyperon is to be produced it must be locally accompanied by an anti-Omega or by an anti-Lambda plus two $`K^+`$ and a $`\pi ^{}`$ (charge, strangeness and baryon number conservation), etc.. Of course this requirement is relatively hard to meet locally, and it inflicts a high statistical ”punishment factor”. We call this effect canonical suppression. In the grand canonical case nothing couples to the Omega locally, we only require three strange quarks and some energy density. The resulting enhancement of the production rate can be understood as an absence of canonical suppression. This statement forms the essence of the analysis  of the transition from C to GC behaviour that we showed in Fig.7.

Before looking at the formal expression for this transition, as derived by Tounsi and Redlich  let me try a crudely oversimplified qualitative argument for illustration, involving the punishment factor only. Consider Omega formation by decay of a partonic medium. Canonical decay then requires 3 s quarks and 3 $`\overline{s}`$ quarks jointly participating (plus other factors that I ignore here). The penalty factor arising from local strangeness and baryon number conservation will, in this naive picture, be $`exp(6m_s/T)`$ whereas in a GC situation it diminishes to $`exp(3m_s/T)`$, the conservation being taken care of elsewhere in the large volume. The enhancement factor then is $`exp(3m_s/T)`$. Taking $`m_s`$=140 MeV and $`T`$=165 MeV we get a factor of about 13: rather close to the ratio of about 15 for the Omega production per participant nucleon pair in central Pb+Pb collisions at top SPS energy, relative to the p+p production rate at similar energy . The Omega is, of course, a very rare species at all incident energies. The ”bulk” strangeness abundance in the system is more directly reflected in the singly strange $`K^+`$ yield. Let us try: the canonical penalty factor would naively be $`exp(2m_s/T)`$ for a $`K^+`$, but only $`exp(m_s/T)`$ in a GC scenario. The GC enhancement factor of bulk strangeness would thus be $`exp(m_s/T)=2.3`$: close enough to the experimentally observed  kaon enhancement, by a factor of about two.

Let us take the above as an amusing illustration concerning the Boltzmann phase space factor for different hadronic species occuring in the denominator of equ. (1). We have, of course, ignored the effect of the baryochemical potential $`\mu _B`$; and many other appropriate considerations. An exact treatment can be found in ref. : Tounsi and Redlich show that the suppression factor governing the yields of strangeness 1, 2, 3 hadrons in proceeding from the canonical case (elementary hadron collisions and $`e^+e^{}`$ annihilation to hadrons) to the GC limit (central AA collisions) is given by

$$F=I_s(V)/I_0(V)$$
(5)
where the $`I`$ are modified Bessel functions of ascending order in strangeness $`S=1,2,3`$ and $`V`$ is the quantum number coherent fireball volume that feeds downward to the frozen-out hadronic gas volume (which may initially be of equal size to the decaying fireball volume). The Bessel functions essentially depend on the coherent volume; the asymptotic limits for this ratio are 1 for $`V\mathrm{}`$, and $`V/2`$ for $`V`$ near zero. The ascend to the ”no suppression limit” at large volume is steepest for singly strange hadrons, and slower for multiply strange species. At small volumes, the factor grows about linearly with the reaction volume. Fig.7 refers to the top SPS energy, $`E/A`$=158 GeV corresponding to $`\sqrt{s}`$=17.3 GeV per participant nucleon pair. From this calculation we infer an enhancement factor of two for singly strange hadrons in the large coherent volume limit, and a much larger factor for $`s`$=2 and 3 hyperons. This study also shows that the canonical suppression - or grand canonical enhancement factor depends, in more detail, on $`\sqrt{s}`$: it falls down with $`\sqrt{s}`$ increasing from SPS to RHIC energies. As we have seen this analysis captures the observations made at top SPS energy . However it also holds at lower SPS energy  and further downward to AGS  and even SIS  energy. Strangeness enhancement (GC behaviour) is, thus, not indicating deconfinement, as proposed by its pioneers . There must be a more general property of the maximum density fireball stage causing GC behaviour: probably large size and coherence. On the other hand let us avoid a misunderstanding: a QCD quark-gluon plasma represents one such case of large size and quantum coherence, and if it decouples directly to a hadron gas one expects to observe GC order . But there ought to be other forms of fireball matter that decouple similarly.

### 5.1 Fireball Size and Coherence

In order to pin down the origin of GC freeze-out let us first look at the size dependence. At top SPS energy we have data for central C+C, Si+Si, S+S and Pb+Pb collisions from NA35 and NA49 . Fig.9 shows compilations of these data  concerning the average total $`\pi ^+,K^+,K^{}`$ and $`\mathrm{\Phi }`$ multiplicities of such central collisions (for Pb+Pb we include the top two centrality windows). The ratio of kaon and $`\mathrm{\Phi }`$ multiplicities to those observed for pions is plotted against the projectile nucleon participant number. We see that all these ratios shoot up steeply from the (canonical) $`pp`$ values that are included for reference . The overall trend resembles the Tounsi-Redlich model  illustrated in Fig 7. We see that S+S already exhausts the GC limit to about 80%, with $`<N_{part}>55`$. In general the model predicts a transition from C to GC order that is yet steeper than the data. However we have to be precise here. Firstly, the model is run with constant temperature and baryochemical potential - perhaps not exactly true in the data. Second, $`<N_{part}>`$ has a different meaning in the data and the model: the statistical model has in reality no such $`<N_{part}>`$ parameter whatsoever, but only a reaction volume in which to enforce either local quantum number conservation or conservation on the average only (from the C to the GC situation). I.e. the authors of Fig.7  fit a coherence volume V, starting from $`pp`$ analysis where they find $`V7fm^3`$ in a situation with, undoubtedly, two participants. With this normalization they merely relabel their volume by $`N_{part}=2V/7fm^3`$. Comparison with AA data is thus not straight forward. We may read off Fig.7 that for s=1 the grand canonical limit is reached for $`\mathrm{"}N_{part}\mathrm{"}20`$ which thus really means $`V70fm^3`$. Not unreasonable: the rms radius of $`{}_{}{}^{32}S`$ is about $`3fm`$, thus its area 19$`fm^2`$. The rapidity distribution of $`K^+`$ has a FWHM of about 2.5 units . Employing the Bjorken-equivalence  which equates the unit rapidity width with a spatial cylinder length of $`1fm`$ we estimate the primordial source volume of a central S+S collision to amount to about $`50fm^3`$.

As a last, important step in our attempt to understand the onset and origin of grand canonical behaviour we shall demonstrate now that it is neither the participant number nor the mere volume of the fireball by itself that matters (recall Fig.2). To this end we show in Fig.10 (left column) the data of Fig.9 confronted with the results of minimum bias PbPb collisions at similar energy. They exhibit a strikingly different, smooth ascend with $`N_{part}`$ which thus can not be the right scaling with variable.

We are now prepared to formulate the essential hypothesis that we have repeatedly hinted at in sections 2-4. We propose that the crucial GC coherence effect resides in the space-time density of ”successive” collisions at the microscopic level, during fireball formation time. This dynamical quantity is, of course, not an ingredient of the statistical model description of hadronization. But it might be the property that justifies a grand canonical description of the system right after hadronization. In order to obtain an estimate of the space-time density of ”successive” collisions in the early fireball, we have to consult a model of the microscopic dynamics. The calculation employs the Frankfurt UrQMD model .

These considerations are illustrated  in the right column of Fig.10. The data of Fig.10 (left column) are replotted here as a function of the average density of successive inelastic collisions that would characterize each of the collisions systems if a quasi classical cascade picture would be applicable. Important: we shall conclude shortly that this picture of sequentiality should not be applicable but we consider the implicit transport ansatz nevertheless to be reliable insofar as it captures the consequences of initial impact geometry and density distributions in a symbolic quantity ”inelastic collisions per cubic fermi and unit time fm/c”, in the emerging fireball. Fig.10 shows a smooth increase of the strangeness enhancement quantities, with saturation occuring in central Pb+Pb. All colliding systems line up smoothly, and we see that saturation, i.e. the GC limit, is reached at $`<\varrho >`$ 6.

Now let us, finally, turn the argument around. A system with 6 separate inelastic collisions concurrent, on average, in each of its space-time unit cells is a strange idea! The uncertainty relation tells us that the energy uncertainty at each ”successive” step would exceed 1.2 GeV here. I.e. the average density of energy uncertainty would exceed the critical QCD energy density. Such a system can not possibly be described as a sequence of hadron collisions. Its decay must occur under global quantum mechanical coherence, from interfering local subprocesses (for which we have no microscopic model as of yet). This mechanism should be the origin of the apparent grand canonical coherence. The cascade-transport model itself thus indicates the point at which it becomes invalid (indeed it does not predict the GC multiplicity pattern). The new physics seems to set in with C+C already, and is well established from S+S onward (Fig.9) where the hot core of the fireball may have a volume of about 50 fm3. Not terribly far from the prediction  of Fig.7 for strangeness one.

It would be extremely interesting to check also the ”onset curves” for s=2,3 with hyperon data for lighter collision systems. Meanwhile we conclude that the occurence of a global chemical potential in the grand canonical description reflects quantum mechanical coherence over a significant volume at the stage of hadronization decay of the fireball. It is important to repeat that in this view strangeness enhancement is not, by itself, an unambiguous ”smoking gun” signal of deconfinement and quark gluon plasma formation. High density hadronic matter (whatever that might be) should also exhibit grand canonical hadron yields. In both cases, however, strangeness enhancement signals the production of hitherto unknown QCD states of matter. This should even apply to the relatively modest effects observed in $`pA`$ collisions  that we mentioned in sect. 2. The ”tube” of target matter struck by the incident nucleon should also be a coherent (albeit small) system.


format: text

text:
# Mode-coupling theory for reaction dynamics in liquids

## I Introduction

Chemical reactions in condensed phases are frequently described using two different methods.Hanggi90  The simple approach consists of a Brownian particle moving in a one-dimensional bistable potential. In this case the dynamics of this system can be described by the Langevin equation.Kramers40  Grote and Hynes Grote80  have extended this approach to non-Markovian processes for a parabolic barrier, where the dynamics are given by the generalized Langevin equation (GLE).Zwanzig61  These and other theories were tested numerically by Straub, Borkovec and Berne,Berne85a ; Berne86  and their results were crucial in the development of the Pollak, Grabert, Hänggi turnover theory.Pollak89

The other approach, which became more common, is based on molecular dynamics simulations, and can account for a general form of the Hamiltonian. Special methods have been developed to accelerate the barrier crossing and thereby numerically determine the rate constant. Chandler Chandler78  showed that in the time-correlation approach to rate constants,Anderson73 ; Bennett77  the reactive flux rapidly decays to a plateau value, which can then be associated with the slow rate for crossing the barrier.Montgomery79  For cases in which no preconceived notion of mechanism or transition state is known, a more general computational method, called the “transition path sampling,” has recently been developed and applied to complex reactive systems.Chandler98 ; Chandler01 ; Chandler02 ; Dellago02

One of the major limitations of the simple approach based on the GLE formalism is related to how one can evaluate the frequency dependent memory friction (related to the random force by the fluctuation dissipation theorem). Most applications assume a Gaussian random force for the GLE. In this limit, the dynamics can be transformed into a Hamiltonian description where the system is linearly coupled to a harmonic bath.Zwanzig73 ; Lindenberg85  Bagchi and his coworkers have used a mode-coupling theory BalucaniZoppi ; BoonYip ; HansenMcDonald  to generalize this approach in order to include non-Gaussian fluctuations.Bagchi97 ; Bagchi99  They assumed that the translational friction is the principal quantity that regulates the diffusive Brownian motion of the reactive system near the barrier region. In this limit, the memory friction is simply given by the Sjögren and Sjölander mode-coupling expression for the memory kernel of a GLE for self-motion in neat liquids.Sjogren79

In this paper we develop an alternative theory to treat non-Gaussian fluctuations within the GLE formalism. First, following the work of Oppenheim and his coworkers,Oppenheim70 ; Oppenheim82 ; Oppenheim91  we derive a GLE for the dynamics of a reactive system in a liquid host. Based on the formal expression for the memory friction of this GLE, and using a combination of kinetic and mode-coupling theories, we obtain a simple expression for the memory friction. This memory friction is then used to obtain the reaction rate within the framework of the Grote-Hynes theory. Our approach is different from that of Bagchi and coworkers Bagchi97 ; Bagchi99  in two ways. First, we need not assume that the translational friction is the principal quantity that regulates the diffusive Brownian motion of the reactive system. Second, it can be shown that this assumption limits the form of the couplings between the system and the liquid host. Thus, a more general form of the Hamiltonian can be treated within our formulation.

Our paper is organized as follows: In Section II we provide a derivation of the GLE, an overview of the Grote-Hynes theory, and a derivation of our kinetic and mode-coupling theories. Tests of our approach for a model isomerization reaction in a LJ fluid is provided in Section III. Finally, we conclude in Section IV.

## II Theory

The development of our approach to reaction dynamics is based on three steps. First, following the work of Oppenheim and his coworkers,Oppenheim70 ; Oppenheim82 ; Oppenheim91  we derive a GLE for the dynamics of a subsystem. Second, we adopt the Grote-Hynes theory,Grote80  which relates the reaction rate to the memory kernel of the GLE. Finally, we develop a theory based on a combination of kinetic and mode-coupling approaches BalucaniZoppi ; BoonYip ; HansenMcDonald  to calculate the reaction rate within the Grote-Hynes formalism.

### II.1 Generalized Langevin Equation

Consider a general Hamiltonian for a reactive system in a liquid host of the from

$$H=H_s(R,P)+H_b(𝐫,𝐩)+\varphi (R,𝐫),$$
(1)
where $`H_s(R,P)`$ is the Hamiltonian of a reactive one dimensional system with phase-space coordinates $`R`$ and $`P`$, and reduced mass $`\mu `$, $`H_b(𝐫,𝐩)`$ is the Hamiltonian of $`N`$ interacting liquid particles with phase-space coordinates $`𝐫`$ and $`𝐩`$, and mass $`m`$, and $`\varphi (R,𝐫)`$ is the coupling between the system and the solvent.

Using the projection operator formalism of Zwanzig and Mori,Zwanzig61 ; Mori65a ; Mori65b  and following the work of Mazur and Oppenheim,Oppenheim70  we now derive reduced equations of motion for the system. First, we define a projection operator $`𝒫`$ that projects out the bath variables:

$$𝒫B=𝑑𝐫𝑑𝐩\overline{\rho }(R,𝐫,𝐩)BB,$$
(2)
where

$$\overline{\rho }(R,𝐫,𝐩)=\rho _b(𝐫,𝐩)\mathrm{exp}\{\beta [\varphi (R,𝐫)w(R)]\}.$$
(3)
In the above equation $`\beta =1/k_BT`$ is the inverse temperature, $`\rho _b(𝐫,𝐩)=\mathrm{exp}\{\beta H_b(𝐫,𝐩)\}/Z_b`$ is the equilibrium distribution function for the isolated bath, $`Z_b`$ is the partition function of the isolated bath, and $`w(R)`$ is the potential of mean force, given by:

$$\mathrm{exp}\{\beta w(R)\}=𝑑𝐫𝑑𝐩\rho _b(𝐫,𝐩)\mathrm{exp}\{\beta \varphi (R,𝐫)\}.$$
(4)
Using the above projection operator it is straightforward to show that the reduced equations of motion for the system are given by the GLE from:Oppenheim70

$$\dot{R}=\frac{P}{\mu },$$
(5)
and

$$\begin{array}{cc}\hfill \dot{P}=\text{e}^{iLt}& \left[\frac{H_s}{R}\frac{w}{R}\right]+K(t)\hfill \\ & \frac{\beta }{\mu }_0^t𝑑\tau \text{e}^{iL\tau }PKK(t\tau ),\hfill \end{array}$$
(6)
where, as noted above, $`\mu `$ is the reduced mass of the system, the $`\mathrm{}`$ is defined in Eq. (2), $`\text{e}^{iLt}`$ is the classical propagator with $`iL=\{H,\mathrm{}\}`$, and the random force $`K(t)`$ is given by

$$K(t)=\text{e}^{i(1𝒫)Lt}\frac{}{R}(\varphi (R,𝐫)w(R)).$$
(7)
Note that the random force involves projected propagation where the classical propagator $`\text{e}^{iLt}`$ is replaced with the projected propagator $`\text{e}^{i(1𝒫)Lt}`$.

The memory kernel of the above GLE is related to the random force by the fluctuation dissipation theorem, and is given by

$$\zeta (t)=\frac{\beta }{\mu }KK(t).$$
(8)
This correlation function is, in general, a function of time and a function of the system coordinate $`R`$ (cf. Eq. (7)).

### II.2 Grote Hynes Theory

To obtain the reaction rate one needs to solve the GLE given by Eq. (6). However, as a result of the fact that the random force involves projected propagation of the couplings between the system and the bath, a complete solution of this GLE is an extremely difficult task. To circumvent this problem, Grote and Hynes have developed an approximate theory to obtain the rate within the GLE formulation.Grote80

The basic assumption made by Grote and Hynes is that in the barrier region, the dynamics can also be described by the above GLE, where the deterministic force (give by $`\left[\frac{H_s}{R}\frac{w}{R}\right]`$) is replaced with an inverted parabolic approximation. After some lengthy algebra, Grote and Hynes obtained the following simple result for the transmission coefficient $`\kappa `$:

$$\kappa =\frac{k}{k_{TST}}=\frac{\lambda _r}{\omega _b}.$$
(9)
The transmission coefficient is given by the ratio of the reactive frequency ($`\lambda _r`$) to the barrier frequency ($`\omega _b`$). The former is related to the memory friction appearing in the GLE, and is given by

$$\lambda _r=\frac{\omega _b^2}{\lambda _r+\widehat{\zeta }^{}(\lambda _r)},$$
(10)
where $`\widehat{\zeta }^{}(\lambda _r)`$ is the Laplace transform of $`\zeta ^{}(t)`$ given by

$$\widehat{\zeta }^{}(\lambda _r)=_0^{\mathrm{}}𝑑t\text{e}^{\lambda _rt}\zeta ^{}(t),$$
(11)
and the symbol “$``$” indicates that the observable is evaluated at the saddle point, namely, in the barrier region as required by the Grote-Hynes theory. To obtain the rate Eq. (10) must be solved self-consistently.

The Grote-Hynes expression for the rate reduces to Kramer’s result if the memory friction decays rapidly on timescales faster than the system’s motion.Grote80  In the strong friction limit, $`\widehat{\zeta }^{}(\lambda _r)/\mu >>\lambda _r`$, the reactive frequency will be much smaller than the barrier frequency resulting in a small transmission coefficient. In the weak friction limit, $`\widehat{\zeta }^{}(\lambda _r)/\mu <<\lambda _r`$, the reactive frequency simply equals the barrier frequency ($`\lambda _r=\omega _b`$), and the rate is given by the transition state rate. This is a shortcoming of the approximation made by Grote and Hynes which fails to capture the turn over to Kramer’s spatial diffusion regime. Hence, the Grote-Hynes rate expression is valid only for high frictions.

The approach developed by Grote and Hynes reduces the complexity of solving the GLE to that of estimating the memory friction of the GLE. For such problems, a very powerful approach has been developed and applied to many interesting dynamical problems in solutions. In the following subsection we describe a combination of kinetic and mode-coupling theories that we have developed to model the memory friction of the GLE in order to study reaction dynamics in liquids.

### II.3 Kinetic and Mode-Coupling Theories

One major difficulty in evaluating the memory kernel of the GLE given in Eq. (6) is that the correlation of the random force involves propagation of $`\varphi (R,𝐫)w(R)`$ in the projected subspace spanned by $`Q=1𝒫`$ (see Eq. (7)). To overcome this difficulty we develop a theory to calculate $`\zeta (t)`$ that is based on a combination of kinetic and mode-coupling theories.BalucaniZoppi  This combination has been used to study density and current fluctuations,BalucaniZoppi ; BoonYip ; HansenMcDonald  solvation and relaxation dynamics,Bagchi99 ; Reichman02b ; Egorov03  and nonlinear spectroscopy Denny01 ; Denny02a ; Denny02b  in classical liquids. A quantum mechanical generalization of this approach has recently been developed and applied successfully to study the dynamic response in quantum liquids.Rabani02a ; Reichman01a ; Rabani02b ; Reichman02a ; Rabani02c ; Rabani02d ; Rabani04

The Grote-Hynes approach requires that the memory friction be approximated in the barrier region. We therefore replace the full memory friction of the GLE ($`\zeta (t)`$) with that approximated at the saddle point ($`\zeta ^{}(t)`$), where (as before) the symbol “$``$” indicates that the position variable of the system is taken at the saddle point. Next, we replace $`\zeta ^{}(t)`$ with an approximate form given by Gotze92

$$\zeta ^{}(t)=\zeta _B^{}(t)+\zeta _{MC}^{}(t),$$
(12)
where $`\zeta _B^{}(t)`$ and $`\zeta _{MC}^{}(t)`$ are the “binary” and “mode-coupling” terms of the memory friction, respectively. The fast decaying binary term is determined from a short-time expansion (to second order in time) of the exact memory friction, and is given by

$$\zeta _B^{}(t)=\zeta ^{}(0)\mathrm{exp}((t/\tau )^2),$$
(13)
where the lifetime $`\tau `$ is given by

$$\frac{1}{\tau ^2}=\frac{1}{2}\frac{\ddot{\zeta }^{}(0)}{\zeta ^{}(0)}.$$
(14)
In the above equations $`\zeta ^{}(0)`$ and $`\ddot{\zeta }^{}(0)`$ are the zero and second time moments of the memory friction, and are given by

$$\zeta ^{}(0)=\frac{\beta }{\mu }\left\{\frac{}{R^{}}(\varphi (R^{},𝐫)w(R^{}))\right\}^2,$$
(15)
and

$$\ddot{\zeta }^{}(0)=\frac{1}{m\mu }\left\{\frac{}{R^{}}\frac{}{𝐫}(\varphi (R^{},𝐫)w(R^{}))\right\}^2.$$
(16)
As noted above, $`\mu `$ is the reduced mass of the system, and $`m`$ is the mass of a liquid particle.

The slow decaying mode-coupling portion of the memory kernel, $`\zeta _{MC}^{}(t)`$, must be obtained from a mode-coupling approach. The basic idea behind this approach is that the random force projected correlation function decays at intermediate and long times predominantly into modes that are associated with quasi-conserved dynamical variables. It is reasonable to assume that the decay of the memory kernel at long times will be governed by those modes that have the longest relaxation time. In the present application the slow decay is basically attributed to couplings between wavevector-dependent density modes of the form

$$b_{𝐤,𝐪}=c_𝐤n_𝐪,$$
(17)
where the self-density mode is given by

$$c_𝐤=\text{e}^{i𝐤𝐫_j},$$
(18)
with $`𝐫_j`$ being the coordinate of liquid particle $`j`$, and the density mode is given by

$$n_𝐤=\underset{j=1}{\overset{N}{}}\text{e}^{i𝐤𝐫_j},$$
(19)
for $`N`$ liquid particles.

In practice, the simplest way to extract the dominant slow contribution of the decay of the memory friction is to introduce another projection operator, $`P_2`$ given by

$$P_2=\underset{\mathrm{𝐤𝐪}}{}\frac{b_{𝐤,𝐪}}{NS(q)}b_{𝐤,𝐪},\mathrm{},$$
(20)
that projects any variables on the space spanned by the slow modes $`c_𝐤`$ and $`n_𝐪`$. Then, following the common approximations in which the projected dynamics of the random force is replaced with the dynamics projected onto these slow variables, and replacing four-point density correlations with a product of two-point density correlations,BalucaniZoppi ; BoonYip ; HansenMcDonald  we find that $`\zeta _{MC}^{}(t)`$ is given by

$$\begin{array}{cc}\hfill \zeta _{MC}^{}(t)& =\frac{\beta }{\mu \rho }\frac{dq^3}{(2\pi )^3}\left[|V^{}(𝐪)|^2/S(q)^2\right]\hfill \\ & \left[F_s(q,t)F(q,t)F_{sb}(q,t)F_b(q,t)\right].\hfill \end{array}$$
(21)
In the above equation $`\rho `$ is the liquid number density, $`S(q)`$ is the structure factor of the neat fluid, $`F_s(q,t)`$ and $`F(q,t)`$ are the self-intermediate and intermediate scattering functions of the neat fluid, respectively, and the vertex is given by the static average of the product of the coupling force and the slow modes:

$$V^{}(𝐪)=\frac{}{R^{}}[\varphi (R^{},𝐫)w(R^{})]b_{𝐪,𝐪}.$$
(22)
As before, the symbol “$``$” indicates that the coupling force is evaluated at the saddle point.

The binary self-intermediate and binary intermediate scattering functions are given by

$$F_{sb}(q,t)=\mathrm{exp}\left\{\frac{q^2t^2}{2\beta m}\right\}$$
(23)
and

$$F_b(q,t)=S(q)\mathrm{exp}\left\{\frac{q^2t^2}{2\beta mS(q)}\right\}$$
(24)
The subtraction of the product of these terms in Eq. (21) is done to prevent over-counting the total memory kernel at short times, namely, to ensures that the even time moments of the total memory kernel are exact to forth order in time.

To obtain the memory friction one requires as input the first two moments $`\zeta ^{}(0)`$ and $`\ddot{\zeta }^{}(0)`$, the static structure factor $`S(q)`$ and the vertex $`V^{}(𝐪)`$. These static averages can be obtained from simulations, or from the proper integral equation formulation. In addition, we need the self-intermediate and intermediate scattering functions. These time-dependent correlation functions of the neat fluid can be obtained from simulations, or alternatively from a similar mode-coupling approach, where a GLE for the density fluctuations is solved using a combination of kinetic and mode-coupling theories.BalucaniZoppi ; BoonYip ; HansenMcDonald

## III Results

To assess the accuracy of the proposed theory we have studied the rate of a model for an isomerization reaction of a diatomic molecule in a Lennard-Jones (LJ) fluid. The Hamiltonian of the entire system and bath can be described by Eq. (1). The isomerizing diatomic molecule is made of two atoms with equal mass $`m^{}`$, interacting via a symmetric double-well potential. Without loss of generality, we place the atoms along the $`z`$-axis, where the position vectors are $`𝐑_1=\{0,0,0\}`$ and $`𝐑_2=\{0,0,R\}`$ for atom $`1`$ and $`2`$, respectively. The reaction coordinate is taken as the distance $`R=|𝐑_1𝐑_2|`$ separating the atoms. We allow the atoms to move only along the $`z`$ axis (which is the reaction coordinate). The reactive system Hamiltonian can be described by:

$$H_s(R,P)=\frac{P^2}{2\mu }+V_0\{4[(RR^{})/a]^21\}^2,$$
(25)
where $`\mu =m^{}/2`$ is the reduced mass of the system, $`V_0`$ is energy barrier separating reactants from products, $`R^{}`$ is the location of the transition state, and $`a`$ is the distance between the two minima corresponding to stable reactants and products. For the results shown below we take $`R^{}=2^{1/6}\sigma `$, $`V_0=5k_BT`$, $`T=2.5ϵ`$, $`a=2\sigma /3`$, and a reduced mass $`\mu =1/2`$ for the isomerizing diatomic molecule.

The solvent Hamiltonian is given by the LJ form:

$$H_b(𝐫,𝐩)=\underset{j=1}{\overset{N}{}}\frac{𝐩_j^2}{2m}+\underset{i>j=1}{\overset{N}{}}4ϵ\left[\left(\frac{\sigma }{r_{ij}}\right)^{12}\left(\frac{\sigma }{r_{ij}}\right)^6\right],$$
(26)
where $`𝐫_j`$ is the position vector of liquid particle $`j`$ with momentum $`𝐩_j`$ and mass $`m`$, and $`r_{ij}=|𝐫_j𝐫_i|`$.

For simplicity we take the solute-solvent interaction to be a two-site LJ potential in which each site of the molecule interacts with the solvent atoms through precisely the same LJ potential:

$$\begin{array}{cc}\hfill \varphi (𝐑,𝐫)& =\underset{j=1}{\overset{N}{}}4ϵ\left[\left(\frac{\sigma }{r_{1j}}\right)^{12}\left(\frac{\sigma }{r_{1j}}\right)^6\right]+\hfill \\ & \underset{j=1}{\overset{N}{}}4ϵ\left[\left(\frac{\sigma }{r_{2j}}\right)^{12}\left(\frac{\sigma }{r_{2j}}\right)^6\right].\hfill \end{array}$$
(27)
Here $`r_{ij}=|𝐫_j𝐑_i|`$. In the application of the mode-coupling theory we need to evaluate the derivative of $`\varphi (𝐑,𝐫)`$ with respect to the reaction coordinate $`R`$. This is done using the chain rule, and the results is given by

$$\frac{\varphi (𝐑,𝐫)}{R}=\frac{1}{2}\left\{\frac{\varphi (𝐑,𝐫)}{Z_2}\frac{\varphi (𝐑,𝐫)}{Z_1}\right\},$$
(28)
where $`Z_1`$ and $`Z_2`$ are the $`z`$ components of the position vector of the diatomic molecule $`𝐑_1`$ and $`𝐑_2`$, respectively, and the factor $`\frac{1}{2}`$ comes from the Jacobian. The derivatives of $`\varphi (𝐑,𝐫)`$ along the other directions vanish since the system is frozen along the $`z`$-axis.

To test our approach, we have calculated the static and dynamic input required to obtain the time-dependent friction, using the Monte Carlo and molecular dynamics simulation techniques. While other approaches, such as the integral equation theory and a proper mode-coupling treatment of density fluctuations, can be used to obtain the static and dynamic input, we feel that a fare test of our theory should rely on numerically exact input. Thus, despite the success of the theoretical approach in predicting structural Rabani01b  and dynamical Rabaniun  properties in LJ systems, we limit our study to the more accurate simulation approach.

The static input was obtained using the molecular dynamics (MD) method for $`N=500`$ particles (including the diatomic) and for a set of densities $`\rho =0.5,0.6,0.7,0.8,0.9,1.0`$ in reduced LJ units.AllenTildesley  Approximately $`10^6`$ MD steps were made for each density. Every $`50`$ steps we collected data for $`\zeta ^{}(0)`$, $`\ddot{\zeta }^{}(0)`$, $`S(q)`$ and for the vertex $`V^{}(𝐪)`$. The self-intermediate and intermediate scattering functions were calculated using the molecular dynamics method for identical conditions. The results for each density were averaged over $`10`$ different runs of total run time of $`t=15`$ reduced LJ units.

In Fig. 1 we plot the static structure factor and the vertex squared for all densities studied in this work. The features observed for the static structure factor are well understood, and have been discussed elsewhere.HansenMcDonald  The magnitude of the vertex given by $`V^{}(𝐪)=\frac{}{𝐑^{}}[\varphi (𝐑^{},𝐫)w(𝐑^{})]b_{𝐪,𝐪}`$ determines the contribution of the different density modes to the decay of the memory friction at intermediate and long times. We find that at high liquid densities the major contribution comes from modes with a characteristic wavelength of $`q2\pi /\sigma `$. This wavelength corresponds to the average interparticle separation, and its value slightly increases with decreasing density (note the shift in the position of the maximum of $`S(q)`$ with density). Thus, density modes associated with liquid motion on length scales of the interparticle separation (namely, on length scale of $`\sigma `$) contribute the most to the friction at intermediate and long times.

As expected, the contribution of lower wavelengths to the decay of the memory friction becomes more significant at lower liquid densities. In fact, at the lowest density studied, the contribution of density modes below $`q2\pi /\sigma `$ to the decay of the memory friction at intermediate and long times is more significant than the contribution of modes near the first peak in $`S(q)`$. We attribute this effect to the change in the mechanism of self-diffusion of liquid particles. At high liquid densities, the self-diffusion is dominated by opening of the cage surrounding a system, with a typical length scale of $`\sigma `$ related to the size of the cage. This is not true at low liquid densities, where liquid particles can hop over much larger distances. In addition, we find that motion within the cage become significant at higher liquid densities, as reflected in non vanishing values of the vertex above $`q2\pi /\sigma `$.

Using the static input obtained from the Monte Carlo simulations along with the values of the self-intermediate and intermediate scattering functions obtained from molecular dynamics simulations, we have calculated the total memory friction given by Eqs. (13) and (21). The results are shown in Fig. (2) for all densities studied in this work. As can be seen, the decay of the memory friction is characterized by two time scales. A fast decay dominated by $`\zeta _B^{}(t)`$ followed by a slower decay dominated by $`\zeta _{MC}^{}(t)`$. The contribution of the mode-coupling portion to the total memory friction is significant only at high liquid densities. At the lower densities studied the mode-coupling portion of the memory friction is negligible.

For the specific model studied here where the solute-solvent interactions equal to the solvent-solvent interactions, the memory friction can also be obtained by inverting the GLE for the velocity autocorrelation function of a neat fluid.Berne70a  Straub, Borkovec, and Berne Berne88  have calculated the velocity autocorrelation function for the LJ fluid using the molecular dynamics technique, and obtained the memory friction by inverting the proper GLE. Our results obtained from the mode-coupling theory agree well with their simulated results (not shown). Specifically, we find that the decay of the memory friction at short times is nearly independent on the liquid density, in agreement with the molecular dynamics results.Berne88  Moreover, at high liquid densities, our theory provides semi-quantitative agreement with the molecular dynamics results at all times. For low liquid densities, we observe small deviations from the simulation results at intermediate times. This shortcoming of the mode-coupling approach is expected since the simulated memory friction becomes slightly negative at intermediate times, and the mode-coupling approximation is known to fail under such circumstances. However, since the contribution of the mode-coupling portion to the memory friction is relatively small at these low densities, this has a vanishing effect on the value of the transmission coefficient.

Using the memory friction obtained from the kinetic and mode-coupling theory we have calculated the transmission coefficient by solving Eq. (9) self-consistently. The results for $`\kappa `$ as a function of the integrated friction $`\widehat{\zeta }^{}(0)`$ (cf. Eq. (11)) are shown in Fig. 3. We compare our results to the results obtained by the weak collision theory (based on a simple connection formula) using an exponential friction,Berne88  and to Kramer’s theory.Hanggi90  For the present model the weak collision theory provides a quantitative agreement for the transmission coefficient (within the noise level of the simulations) in comparison with results obtained using the absorbing boundary approximation Berne85b ; Berne85c  to the reactive-flux formalism.Chandler78 ; Berne85  As can be seen, our theory provides quantitative results for the transmission coefficient over the entire range of frictions studied.

The different values of the integrated friction $`\widehat{\zeta }^{}(0)`$ at which we have calculated $`\kappa `$ where obtained by scanning the liquid density from low $`\rho =0.5`$ to high $`\rho =1.0`$ values. For each density we have calculated the static and dynamic input required by our mode-coupling theory, and obtained the transmission coefficient by solving Eq. (9) self-consistently. Alternatively, one can control the integrated friction $`\widehat{\zeta }^{}(0)`$ by changing the value of the reduced mass $`\mu `$ of the isomerizing diatomic molecule (cf. Eq. 8). We find that our mode-coupling theory provides similar quantitative results (not shown) for the transmission coefficient when the reduced mass of the isomerizing diatomic molecule is varied.

## IV Conclusions

We have presented a theoretical approach for the calculation of reaction rates in condensed phases based on the Grote-Hynes formalism. A combination of kinetic and mode-coupling theories were developed to obtain the memory friction required by the Grote-Hynes formalism. The approach was applied to study a model isomerization reaction of a diatomic molecule in a LJ fluid. Good agreement for the transmission coefficient was obtained in comparison with the simulation results of Straub, Borkovec, and Berne Berne88  based on the reactive flux formalism.

Unlike simulation techniques, our approach is a theory and thus provides additional insight into the reaction dynamics in liquids. For example, we showed that the contribution of the mode-coupling portion to the decay of the memory kernel at intermediate times is significant only at high liquid densities. Thus, an accurate description of the reaction dynamics at low liquid densities can be obtained from a kinetic theory alone. Furthermore, the mechanism for the decay of the memory friction (which is reflected in the reaction rate) is quite different at low versus high liquid densities. At high liquid densities the decay of $`\zeta ^{}(t)`$ is dominated by liquid modes with a length scale comparable to the separation between the fluid particles, while at lower liquid densities motion on larger length scales also contributes the decay of $`\zeta ^{}(t)`$. This is significant for the development of coarse grained models for reaction dynamics.

We believe that our approach will be useful in other situations in which simulation techniques are still limited. For example, for reaction dynamics in supercooled liquids that are characterized by slow density fluctuations. Or for liquid hosts that are characterized by quantum mechanical susceptibilities. Work along these directions is currently underway.

## V Acknowledgments

This work was supported by The Israel Science Foundation founded by the Israel Academy of Sciences and Humanities (grant number 31/02-1). The authors would like to thank David Reichman for stimulating discussions and Oded Hod for support with the simulations.