text:
# Ab-initio determination of the localized/delocalized 𝑓-manifold in UPd2Al3

## Abstract

The electronic structure of UPd2Al3 is described using the self-interaction corrected local-spin-density approximation to density functional theory. The groundstate is found to be characterized by the coexistence of localized ($`f^2`$) and delocalized U $`f`$ electrons, in agreement with experimental evidence. We observe significant difference in electronic structure between UPd2Al3 and the previously studied UPt3 compound. Even though a trend towards localization exists in UPt3, the total energies and the density of states at the Fermi level favor a groundstate with localized $`f^1`$, rather than $`f^2`$ U ions.

Heavy fermion (HF) systems are some of the most fascinating condensed matter materials. Their highly correlated electrons give rise to large values of the low temperature specific heat coefficient, $`\gamma `$, which in certain cases is larger than that of a normal metal by more than a factor of 1000. This is associated with highly enhanced effective quasiparticle masses, $`m^{}`$, of some 1000 times the free electron mass. Recently new insight has been gained into the nature of these heavy quasiparticles in the magnetically ordered heavy fermions UPd2Al3 and UPt3. This work has established that magnetic order is the prerequisite to formation of these heavy quasiparticles. Specifically, it is asserted that coupling of localized and delocalized 5f electrons lies at the origin of heavy quasiparticles in UPd2Al3 and UPt3.

In the standard description, the large specific heat coefficients can be correlated to a large density of states (DOS) of heavy quasiparticles at the Fermi level, i.e., a narrow resonance at the Fermi level resulting from the interaction of the localized $`f`$ states with the conduction electrons. In the weak-coupling Kondo limit the atomic moments of the localized states are individually screened by the conduction electrons. With increasing hybridization, coherence sets in, leading to the Fermi liquid of weakly interacting, extremely heavy quasiparticles. However, the heavy fermion state will not form, if the interactions between unscreenend or incompletely screened local moments induce a magnetic phase transition, with the onset of magnetic order. This is the case in many rare earth and actinide metals such as e.g. UPd3. Contrary to this picture, the new observation by Dressel  indicates that for UPd2Al3 and UPt3 the heavy quasiparticles are only created by magnetic excitations in the magnetically ordered state, within the manifold of localised and delocalised f states. In the present paper these localised and delocalised states in UPd2Al3 are determined from ab initio self-interaction corrected (SIC) local spin density (LSD) calculations and discussed in comparison with UPt3.

In UPd2Al3, the specific heat coefficient is $`\gamma `$=140 mJ mol-1 K-2 i.e. less extreme than the value for the prototypical heavy fermion compound CeAl3 ($`\gamma `$=1600 mJ mol-1 K-2) but still considerably larger than the value for a normal metal such as Pd ($`\gamma `$=9.4 mJ mol-1 K-2). At TN=14.3 K, UPd2Al3 undergoes an antiferromagnetic phase transition, with a resulting moment of 0.85 $`\mu _B`$ per U-ion. At the critical temperature Tc=2 K,  a superconducting phase transition occurs, as indicated by a sharp drop in electrical resistivity, and a large jump in the specific heat.  The pronounced magnetic moment is retained below Tc From the specific heat measurements,  muon spin rotation measurements,  and neutron scattering studies,  it emerges that superconductivity and anti-ferromagnetism coexist in UPd2Al3. It has been established that magnetic excitations between two $`f`$-electrons localized at each site, and the remaining $`f`$-electron are causing superconductivity.  Recently it has been suggested  that the heavy fermion superconductor UPt3 can be similarly described in terms of a coupling between the localized $`f^2`$ configuration and the delocalized $`f`$’s.

Applying density functional theory  (DFT) to such HF compounds is not straightforward. While this theory is formally exact for an arbitrary many-electron system, in its most widely used approximations, like the local spin density or generalized gradient approximation (GGA), it relies on the exchange and correlation effects of a homogeneous electron gas, which is a rather poor reference system for heavy fermion compounds. The LSD scheme has been very successful in describing the cohesive properties of solid state systems with itinerant valence electrons. In the 4$`f`$ and 5$`f`$ systems, however, the band formation competes with a strong tendency towards localization due to large on-site $`f`$-$`f`$ correlations that go beyond the LSD. Therefore, it is not obvious to what extent the results from bandstructure calculations remain valid. Nevertheless, owing to the fact that the 5$`f`$ electrons in the U compounds are considerably less localized than for example the 4$`f`$ electrons in the corresponding rare-earth compounds, LSD is known to provide useful electronic structure in the normal phase of these compounds. For instance, the Fermi surface of UPt3 is well reproduced by LSD, and for UPd2Al3 it was shown  that the experimental de Haas van Alphen results, and the measured magnetic moment, could be reproduced by calculations treating all the U $`f`$-electrons as itinerant, hybridizing with the conduction electrons. The latter work  concluded, however, that in order to explain the very anisotropic magnetic susceptibility of UPd2Al3, one needed to consider the coexistence of two different kinds of 5$`f`$ electrons, i.e., localized and itinerant ones.

The self-interaction-corrected local-spin-density approximation (SIC-LSD),  constitutes a scheme capable of treating both localized and delocalized electrons on an equal footing, by including an explicit energy contribution for an electron to localize (the self-interaction correction).  This localization energy is the sum of the self-Coulomb and self-exchange-correlation energies of the f-state. For an extended state this energy correction vanishes, but for a localized U 5$`f`$-state the self-interaction correction is of the order of 30 mRy. In the SIC-LSD application to actinide systems, a particular localized $`f^n`$ configuration is assumed on each actinide ion, while all other degrees of freedom are allowed to form bands, in a self-consistency cycle. The localized states loose their band-formation energy, but gain the self-interaction correction energy which is subtracted from the LSD total energy functional. On balance, comparing different localized $`f^n`$ scenarios, the global energy minimum establishes which configuration is the groundstate. In the SIC-LSD picture, the localized electrons are no longer available as band-forming valence electrons, and accordingly the valency of the actinide ion may be defined as $`N_{val}`$=$`Z`$-$`N_{core}`$-$`N_{semi}`$-$`N_{SIC}`$, where Z is the atomic number (92 for U), $`N_{core}`$ is the number of atomic core electrons (78), $`N_{semi}`$ is the number of semicore electrons (8), and $`N_{SIC}`$ is the number of localized electrons. Thus, for UPt3 the SIC-LSD scheme has revealed a groundstate characterized by localized $`f^1`$ electron on U (pentavalent configuration), which is intermediate between the itinerant ($`f^0`$) behavior in URh3 and the more localized $`f^2`$ U ions found in UPd3 and UAu3 (tetravalent configurations). In the present work a similar study of UPd2Al3 is reported. Three configurations have been investigated, assuming respectively localized $`f^2`$, $`f^1`$ and $`f^0`$ configurations on the U sites.

The SIC-LSD formalism has been implemented  within the tight-binding linear muffin-tin orbitals method (LMTO),  with the atomic sphere approximation (ASA). All relativistic effects, including spin-orbit interaction, have been incorporated in the Hamiltonian. For UPd2Al3 two separate energy panels have been used, one for the U 6$`s`$ and 6$`p`$ semi-core states, and one for the valence electrons. Apart from the 4$`d`$, 5$`s`$ and 5$`p`$ orbitals on the Pd atom, and the 3$`s`$, 3$`p`$, and 3$`d`$ orbitals on the Al atom, the valence panel has comprised the 7$`s`$, 6$`d`$ and 5$`f`$ orbitals on the U atom with the 7$`p`$ degrees of freedom treated by downfolding.

In figure 1 we show the results of the total energy calculations as a function of volume for respectively the hexavalent ($`f^0`$), pentavalent ($`f^1`$), and tetravalent ($`f^2`$) configurations, with the $`f^0`$ simply corresponding to the standard LSD approximation. The tetravalent groundstate configuration is seen to give the lowest total energy, in agreement with the U4+ configuration predicted by experiment. The global minimum occurs at a volume of 681 au3, in good agreement with the measured value, Vexp=700 au3. The minima for the $`f^0`$ and $`f^1`$ configurations are shifted slightly to lower volumes, as the partial pressure of itinerant $`f`$-electrons is negative. The exact values of the minimum volumes depend somewhat on the ratio of atomic spheres adopted in the ASA, but the relative ordering of the minima is preserved, i.e., the finding of the ground state as $`f^2`$ is robust in this respect. The present calculations use the fixed ratio of atomic spheres as $`R_U:R_{Pd}:R_{Al}=3.4:2.4:2.5`$.

With two of the $`f`$-electrons localized in the tetravalent configuration, the DOS of UPd2Al3 around the Fermi level differs considerably from the corresponding DOS in the LSD configuration, as shown in Fig. 2. Thus, in the itinerant scenario, with all $`f`$-electrons delocalized (Fig. 2a), we obtain a relatively large partial $`f`$-DOS (dotted line) at the Fermi level of about 102 states/Ry, as compared to the $`f^2`$ configuration with only 20 states/Ry (Fig. 2b). In the latter configuration, most of the $`f`$ weight has been transferred into the two localized states, and the remaining delocalized $`f`$-states hybridize with the $`s`$, $`p`$ and $`d`$ levels on the neighbouring Al and Pd sites. This strong hybridization between $`f`$-electrons and conduction electrons which is observed in both the tetravalent and delocalized  configurations, can explain the fact that the measured specific heat coefficient is only about 140 mJ mol-1 K-2 in UPd2Al3, i.e., we are dealing with moderately heavy itinerant electrons at the Fermi level. Despite this qualitative agreement, we must however not forget that by allowing an $`f`$-electron to delocalize the corresponding strong on-site $`f`$-$`f`$ correlations are ignored, and that we therefore might considerably overestimate the effect of $`f`$-conduction electron band-formation. The effect of these correlations becomes clear when we compare the calculated and experimental values for the specific heat coefficient. Here, in the tetravalent configuration we obtain $`\gamma `$=12 mJ mol-1 K-2, i.e. too small by a factor of 10 with respect to experiment. The same discrepancy between calculation and measurement has been observed in other heavy fermion compounds,  the reason being that although the correlations are actually not strong enough to localize the $`f`$-electrons, they lead to a considerable narrowing of the hybridized $`fpd`$-band, resulting in a specific heat coefficient that is enhanced with respect to the LSD calculations by a factor of up to 20-30. Thus we do expect to find too small a value for $`\gamma `$, and the fact that LSD value is in better agreement ($`\gamma _{LSD}`$=23 mJ mol-1 K-2) with experiment is a result of the increased $`f`$-weight at the Fermi level due to the delocalization of the $`f`$-electrons, rather than an indication of a better description of the electronic structure in terms of narrow bands.

There are now strong indications that, in UPd2Al3, the existence of two localized $`f`$-electrons is the prerequisite for the appearance of superconductivity and heavy quasiparticles.  The argument is that the itinerant heavy quasiparticles are strongly coupled to the localized $`f`$ electrons, and interact effectively with each other by the intermediary of magnetic exitations. It has been suggested that the superconductivity in UPt3 might be of similar origin. Optical conductivity measurements by these authors show a pronounced pseudogap in the low temperature response, in both UPd2Al3 and UPt3, which they relate to magnetic correlations between localized and delocalized $`f`$-electrons. Zwicknagl et al.  have shown that treating two of the $`f`$-electrons in UPt3 as localized in LDA band-structure calculations, may reproduce the observed de Haas-van Alpen frequencies, and that the interaction between itinerant and localized $`f`$-electrons may lead to mass enhancements of a factor $`10`$. There are however noticeable differences in the low temperature properties of UPd2Al3 and UPt3. Thus even though antiferromagnetism and superconductivity coexist in UPt3 below Tc=0.55 K, the magnetic moment is very small, only 0.02 $`\mu _B`$. Also, unlike in UPd2Al3, in UPt3 the magnetic moment is affected by the superconducting phase transition.

The SIC-LSD calculations of UPt3  also find indications in favour of coexisting localized and delocalized $`f`$-electrons. However in those calculations, it is the pentavalent $`f^1`$ configuration, rather than the tetravalent $`f^2`$ configuration, which turns out to be energetically most favourable. Furthermore, in the pentavalent configuration a narrow $`f`$-peak is pinned to the Fermi level, in agreement with the heavy fermion character of this compound, whilst in the tetravalent configuration the UPt3 DOS at the Fermi energy vanishes, as is shown in Fig. 3. In fact, according to the present SIC-LSD calculations, in the scenario with two localized $`f`$-electrons, UPt3 becomes very similar to UPd3, which is then difficult to reconcile with experimental evidence,  since UPd3 is neither a heavy fermion compound nor a superconductor. In UPd2Al3, as can be seen from Fig. 1, the energy difference between the $`f^2`$ and $`f^1`$ configuration is $`E(f^2)E(f^1)6`$ mRy. For UPt3, UPd3 and UAu3 the corresponding energy differences are $`+5`$ mRy, $`5`$ mRy and $`9`$ mRy, respectively. Thus with respect to localization, UPd2Al3 is quite similar to UPd3 (tetravalent configuration), and more localized than UPt3. Experimental evidence seems to suggest that UPt3 is close to an antiferromagnetic instability,  and consequently the $`f`$-electrons are situated at the borderline of a localization/delocalization transition. In our SIC-LSD calculations for UPt3, this instability is reflected by a groundstate configuration with one extra delocalized $`f`$-electron as compared to UPd2Al3. However, correlations of this extra delocalized electron are most likely not adequately accounted for by the LSD approximation, and would plausibly need to be described in a theory involving fluctuating $`f^1`$ and $`f^2`$ configurations, such as the dynamical coherent potential approximation. In particular, the nearly vanishing magnetic moment observed for UPt3 is not reproduced in the SIC-LSD approach. The SIC-LSD results indicate however that it is energetically not favourable to statically localize two $`f`$ electrons. Hence, the quasiparticle model of Zwicknagl et al.  might possibly still be valid for UPt3, although it seems to be more straightforwardly applicable to UPd2Al3.

In summary, SIC-LSD total energy and DOS calculations of UPd2Al3 have been performed. A ground state of localized U $`f^2`$ shells has been obtained where additionally strongly hybridized states of $`f`$ character appear at the Fermi level. The SIC-LSD ab initio calculations provide direct evidence for the picture of coexisting localized and delocalized $`f`$-electrons which is inferred from experiments (such as ) to explain heavy quasiparticles and superconductivity. The calculations do not confirm a similar degree of localization for the $`f`$-electrons in UPt3 and therefore magnetic excitations in UPd2Al3 and UPt3 are most likely to be different.

This work has been partially funded by the Training and Mobility Network on ‘Electronic Structure Calculation of Materials Properties and Processes for Industry and Basic Sciences’ (contract:FMRX-CT98-0178) and the Research Training Network ’Psi-k f-electron’ (contract:HPRN-CT-2002-00295).


format: text

text:
# The electrical conductivity and soft photon emissivity of the QCD plasma

## Abstract

The electrical conductivity in the hot phase of the QCD plasma is extracted from a quenched lattice measurement of the Euclidean time vector correlator for $`1.5T/T_c3`$. The spectral density in the vicinity of the origin is examined using a method specially adapted to this region, and a peak at small energies is seen. The continuum limit of the electrical conductivity, and the closely related soft photon emissivity of the QCD plasma, are then extracted from a fit to the Fourier transform of the temporal vector correlator.

preprint: TIFR/TH/03-01, hep-lat/0301006

The soft photon production rate from the plasma phase of hadronic matter is of importance to searches for the QCD phase transition wa98 . Consequently, there has been a long history of attempts at perturbative computations of this rate history . The first lattice prediction of dilepton (off-shell photon) rates was performed a while back dilepton . Recently the leading order computation of the photon production rate was completed amy . The Kubo formula relates the soft limit of this rate to the DC electrical conductivity of the QCD plasma, $`\sigma `$. To leading log accuracy in the gauge coupling, $`g=\sqrt{4\pi \alpha _S}`$, one has $`\sigma \alpha T/g^4\mathrm{log}g^1`$, where $`\alpha `$ is the fine structure constant. The proportionality constant has been computed recently in the leading-log approximation amy2 . Here we report the first computation of $`\sigma `$ and the soft photon emissivity from a quenched lattice computation in a region of temperature where $`g`$ is large and the weak-coupling approach fails amynote . Our methods can also be applied to other transport problems.

The photon emissivity at temperature $`T`$ is related to the imaginary part of the retarded photon propagator, i.e., the spectral density, $`\rho _{EM}^{\mu \nu }`$, for the electromagnetic current correlator, through the relation

$$\omega \frac{d\mathrm{\Omega }}{d^3𝐩}=\frac{1}{8\pi ^3}n_B(\omega ,T)\rho _{EM}^{}{}_{\mu }{}^{\mu }(\omega ,𝐩,T).$$
(1)
In this work we shall take $`\omega =𝐩=0`$, and hence obtain the soft photon production rate. Since the EM Ward identity gives $`\rho ^{00}(\omega ,\mathrm{𝟎},T)=0`$, this soft limit is related to transport properties of the QCD plasma through the Kubo formula,

$$\sigma (T)=\frac{1}{6}\frac{}{\omega }\rho _{EM}^{}{}_{i}{}^{i}(\omega ,\mathrm{𝟎},T)|_{\omega =0},$$
(2)
where the sum is over spatial polarisations. A lattice determination of this rate proceeds from the spectral representation for Euclidean current correlators—

$$G_{EM}(t,𝐩,T)=_0^{\mathrm{}}\frac{d\omega }{2\pi }K(\omega ,t,T)\rho _{EM}(\omega ,𝐩,T),$$
(3)
where the integral kernel $`K=\mathrm{exp}(\omega t)n_B(\omega ,T)+\mathrm{exp}(\omega t)[1+n_B(\omega ,T)]`$. $`G_{EM}`$ is the product of the vector correlator summed over all polarisations, $`G_V`$, and the EM vertex factor $`C_{EM}=4\pi \alpha _fe_f^2`$, where $`e_f`$ is the charge of a quark of flavour $`f`$. On discretising the integral it becomes clear that the extraction of $`\rho _{EM}`$ from the lattice computation of $`G_{EM}`$ is akin to a linear least squares problem. The complication is that the (potentially infinite) number of parameters to be fitted exceeds the number of data points (which is half the number of lattice sites in the time direction, $`N_t`$). The solution is to constrain the function $`\rho _{EM}`$ through an informed guess baym , and use a Bayesian method to extract it nullsp . The Maximum Entropy Method (MEM) gubernatis ; asakawa  along with a free-field theory model of the spectral function has been used in the past dilepton . The hard dilepton rate for $`\omega /T4`$ is fully under control, with lattice and perturbation theory in good agreement dilepton . For that reason we concentrate here on the electrical conductivity and the soft photon rate.

Correlators were investigated at $`T=1.5T_c`$, $`2T_c`$ and $`3T_c`$ in quenched QCD. The temperature range is realistic for heavy-ion collisions. However, $`g>1`$ in this entire range of temperature precise  and is therefore ineffective in the separation of length scales upon which weak-coupling approaches depend. In order to make continuum extrapolations, the computations were performed on a sequence of lattice spacings, $`a=1/8T`$, $`1/10T`$, $`1/12T`$ and $`1/14T`$, (i.e., $`N_t=8`$, 10, 12 and 14). Quark mass effects were controlled by working with staggered quarks of masses $`m/T_c=0.03`$ and $`0.1`$. Details of the runs, statistics, and the generation of configurations for $`N_t<14`$ are described in gavai . For these lattice spacings the computations were performed on two different spatial volumes in order to control finite volume effects. For $`N_t=14`$ we have added runs on $`14\times 30^3`$ lattices for $`T=1.5T_c`$ and $`2T_c`$, and on $`14\times 44^3`$ lattices for $`T=3T_c`$, generating 50 configurations separated by 500 sweeps each. We have measured vector correlators with two degenerate flavours of quarks. It has been demonstrated recently that in this limit the charged and uncharged vector correlators are identical nf11 .

Small but statistically significant differences between the lattice results and ideal quark gas predictions for $`G_{EM}`$ are observed at all temperatures, lattice spacings, quark masses and volumes investigated. In any lattice computation, we expect the high frequency part of $`\rho _{EM}`$ to contain lattice artifacts. Moreover, physics at momenta of order $`1/a`$ is perturbative short  and not of interest in the present context. We remove this physics by taking the difference between the Euclidean temporal propagators in QCD and an ideal quark gas (free field theory) on the same lattice—

$$\mathrm{\Delta }G_{EM}(\omega )=G_{EM}^{QCD}(\omega )G_{EM}^{ideal}(\omega ),$$
(4)
A subtraction is needed to remove the $`\omega ^2`$ divergent pieces from the dispersion relations hilg . Other benefits accruing from this are discussed later.

We have estimated the spectral density by two classes of methods. The first class of general techniques consist of discretising the integral in eq. (3) into $`N_\omega `$ energy bins and rewriting the equation in the form $`G_{EM}=K\rho _{EM}`$ where $`K`$ is now an $`N_t\times N_\omega `$ matrix, $`G_{EM}`$ the data vector of length $`N_t`$ and $`\rho _{EM}`$ represents a vector of length $`N_\omega `$ integ . For $`N_\omega >N_t`$ the solution is non-unique. Additional constraints, called priors, must then be imposed to determine them ill . The extraction of the spectral density is performed in the context of Bayesian parameter extraction. Given the data on $`G`$, the probability distribution function for $`\rho `$ can be written using Bayes’ formula

$$P(\rho |G)=P(G|\rho )P(\rho )/P(G),$$
(5)
where $`P(A|B)`$ denotes the conditional probability of $`A`$ given $`B`$. The probability distribution function $`P(\rho )`$ contains the prior information on $`\rho `$ that is needed for the analysis. Writing $`P(\rho |G)=\mathrm{exp}[F(\rho )]`$, the maximum likelihood analysis of the probability reduces to the problem of minimizing $`F`$. Since $`P(G)`$ is independent of $`\rho `$, this problem can be formulated as a minimisation of the function

$$F(\rho )=(GK\rho )^T\mathrm{\Sigma }^1(GK\rho )+\beta U(\rho ),$$
(6)
where the first term is the logarithm of $`P(G|\rho )`$ and $`P(\rho )=\mathrm{exp}[\beta U]`$. The superscript $`T`$ denotes a transpose, $`\mathrm{\Sigma }`$ is the covariance matrix of the data, and $`\beta `$ is a non-negative parameter whose choice is specified later.

The MEM technique consists of choosing some vector $`\rho _{EM}^0`$ and defining $`U(\rho _{EM})=_i\rho _{EM}^i[\mathrm{log}(\rho _{EM}^i/\rho _{EM}^{0i})1]`$, where the sum is over components of the vectors. In previous works the prior $`\rho _{EM}^0`$ has been chosen to be the vector spectral function in an ideal quark gas dilepton . Another whole class of techniques is obtained by choosing $`U(\rho _{EM})=[(\rho _{EM}\rho _{EM}^0)]^2`$ where $``$ is a non-singular matrix. The choices $`=1`$, $`D`$ and $`D^2`$ (where $`D`$ is a discretisation of the derivative) have been suggested in the literature. $`=1`$ is the model that $`\mathrm{\Delta }\rho _{EM}\rho _{EM}\rho _{EM}^0=0`$ except as forced by the data, $`=D`$ makes the a priori choice that $`\mathrm{\Delta }\rho _{EM}`$ is constant and $`=D^2`$ is the prior choice of smooth $`\mathrm{\Delta }\rho _{EM}`$ honerkamp .

Such regulators have the added advantage that minimisation of the function $`F`$ in eq. (6) yields the linear problem—

$$\left[K^T\mathrm{\Sigma }^1K+\beta ^T\right]\rho _{EM}=K^T\mathrm{\Sigma }^1G_{EM}+\beta ^T\rho _{EM}^0.$$
(7)
Since $`^T`$ is positive definite, it is clear that the term in $`\beta `$ on the left hand side regulates the problem, by adding a term to $`K^T\mathrm{\Sigma }^1K`$ which makes the sum invertible. Since, for a well-determined parameter fitting problem, the value of $`F`$ is the $`\chi ^2`$ value, we choose a value of $`\beta `$ at which $`F=N_t`$ at the minimum of $`F(\beta )`$, i.e., at the maximum a posteriori probability.

Considering the Bayesian problem as a field theory for the function $`\rho _{EM}`$, the method of maximising the a posteriori probability is equivalent to a semi-classical solution. The advantage of choosing a linear regulator is two fold. First, the search for the minimum is simply the solution of a system of linear equations; in non-linear minimisation it is no simple matter to correctly identify the global minimum press . Second, the linear problem is guaranteed to have a single minimum, whereas a general non-linear regulator may have multiple local minima, leading to complications analogous to the physics of phase transitions.

Since previous work has demonstrated that for $`\omega T`$ lattice computations match perturbation theory dilepton , we focus our attention on the region $`\omega \pi T`$. The linear relation between $`G_{EM}`$ and $`\rho _{EM}`$ means that we can assume $`\rho _{EM}=\rho _{EM}^0+\mathrm{\Delta }\rho _{EM}`$, where $`\rho _{EM}^0`$ is the usual MEM prior in an ideal quark gas. At small $`\omega `$ this goes to zero faster than linearly in $`\omega `$ and hence does not contribute to $`\sigma /T`$ deltafoot . By choosing to work with $`\mathrm{\Delta }G_{EM}`$, this $`\rho _{EM}^0`$ is removed from the problem, and we are freed to concentrate on the piece $`\mathrm{\Delta }\rho _{EM}`$, which contains all the information needed to extract $`\sigma `$. Then in eq. (7) we use $`=1`$, replace $`G_{EM}`$ by $`\mathrm{\Delta }G_{EM}`$ and $`\rho _{EM}`$ by $`\mathrm{\Delta }\rho _{EM}`$, and remove the term in $`\rho _{EM}^0`$. The upper limit of the integral was truncated to $`\omega =2n\pi T`$ and the range divided into an uniform mesh of $`N_\omega `$ points. Varying $`n`$ and $`N_\omega `$ independently in the range $`2n4`$ and $`16N_\omega 64`$ has no effect on the quality of the fit to the data (see Figure 1).

Statistical errors on $`\mathrm{\Delta }\rho _{EM}`$ are assigned by a bootstrap over the measured values of $`\mathrm{\Delta }G_{EM}`$. These are minor compared to uncertainties arising from algorithmic parameters. The latter are estimated by changing the integration method which is used to discretise eq. (3), the bin size in the integration, the integration limit, $`N_\omega `$, and the Bayesian prior specified by the matrix $``$. Thrice the combined uncertainty due to these four factors is shown as the band in Figure 2 within which $`\mathrm{\Delta }\rho _{EM}`$ lies. Even with this generous allowance for uncertainties there seems to be a peak in the spectral function at small $`\omega `$. For $`\omega \pi T`$ the spectral function is roughly consistent with free field theory, but there is some evidence of a further peak at $`\omega 16\pi T`$. The most important systematic uncertainty turns out to be related to control over the limit $`a0`$. We found that the position of the peak and the slope at the origin change in going from $`N_t=8`$ to 14. This phenomenon has been noticed earlier in the context of MEM limit . A method which allows for better control of the continuum limit is required.

For this we utilize a second class of Bayesian methods, in which the prior is a model of the observed bump in the soft part of the spectrum. Since $`\rho _{EM}`$ is real and non-singular for real $`\omega `$, odd in $`\omega `$, and non-negative for $`\omega >0`$, one can choose to work with the most general form which gives rise to a non-vanishing electrical conductivity,

$$\frac{1}{T^2}\mathrm{\Delta }\rho _{EM}(z)=\frac{z_{n=0}^N\gamma _nz^{2n}}{1+_{m=1}^M\delta _mz^{2m}},$$
(8)
where $`z=\omega /T`$ and with all $`\gamma _n`$ and $`\delta _m`$ real and non-negative aarts . The constraint that $`\mathrm{\Delta }\rho _{EM}0`$ at large $`\omega `$ is imposed by choosing $`M>N`$ param . We shall use the notation $`(N,M)`$ to denote a particular choice of $`N`$ and $`M`$. Bayesian techniques for parameter estimation then proceed by choosing a priori probability distributions for each parameter lepage .

The parameters in eq. (8) are most conveniently extracted by fitting to the Fourier transform of eq. (3) over the Euclidean time $`t`$ dolan

$$𝒢_{EM}(\omega _n,𝐩,T)=\frac{d\omega }{2i\pi }\frac{\rho _{EM}(\omega ,𝐩,T)}{\omega \omega _n},$$
(9)
where the Euclidean frequencies are $`\omega _n=2in\pi T`$, ($`1nN_t`$ on a lattice) and the path of integration over complex $`\omega `$ runs over the real line and is closed in the upper half-plane. The form of $`\rho _{EM}`$ in eq. (8) can then be used to express the Fourier coefficients in terms of the parameters, which can be determined either through a least squares method if $`1+N+M<N_t/2`$, or a Bayesian method when there are more parameters than data.

A particular simplification occurs for $`(0,M)`$, since $`\gamma _0`$ is the only parameter that contributes for $`\omega _n=0`$. In all these cases $`\sigma /T=C_{EM}(\chi _V\chi _V^0)/3T^2`$, where $`\chi _V=𝒢_V(0,\mathrm{𝟎},T)`$ is the vector meson susceptibility gupta  obtained in QCD and $`\chi _V^0`$ is the same quantity for an ideal quark gas on the same lattice arvind . Since the remaining parameters do not appear in this expression, their prior probabilities can be integrated out of the problem, without any assumptions about them. Such a marginalisation of the prior distribution is a general technique which has been demonstrated on other problems in the past for .

The extraction using $`N=0`$ must be insufficient for $`T>10^{11}T_c`$, since it does not reproduce the parametric dependence of $`\sigma `$ on $`g`$ in weak-coupling theory, which is expected to work at these temperatures. This can be improved by allowing for other values of $`(N,M)`$. We have investigated the stability of our results by going to $`(1,2)`$. Such a multiparameter fit moves the result up by 7%, which is within the statistical uncertainty. The electrical conductivity is thus reasonably stable, although it would be interesting in future to investigate its stability further, especially by using larger values of $`N_t`$.

In principle, such an extraction of parameters other than $`\gamma _0`$ in eq. (8), allows us to proceed beyond the $`\omega =0`$ limit of the dilepton rate. As more parameters are determined, the shape of the soft dilepton spectrum is also better constrained. An interesting open question is of the number of Fourier coefficients needed to fix the shape of the dilepton spectrum. This question is related to the stability of the transport coefficient, and we plan a study in the near future to address this question.

The estimates of $`\sigma `$ from the formula above are subject to lattice artifacts of order $`a^2`$ coming from $`𝒢_V`$. The continuum limit can then be obtained by an extrapolation in $`1/N_t^2`$. Finite volume effects turn out to be invisible within errors. Nor is there any visible quark mass dependence for small quark masses, since $`m/T_c=0.03`$ and $`0.1`$ give identical results within errors. We estimate $`\sigma /T7C_{EM}`$ in the continuum limit of the temperature range we studied. Finally, have used the estimate of $`\sigma `$ to predict the soft photon emissivity of the QCD plasma in equilibrium, as shown in Figure 3.

In summary, we adopted a sequence of Bayesian techniques for the inverse problem of extracting spectral densities, $`\rho _{EM}`$. In view of the results of dilepton  we used the dispersion relations for Euclidean propagators after subtraction of the ideal gas values of the Euclidean temporal correlator, $`\mathrm{\Delta }G_{EM}`$. We observed that in the QCD plasma, in the temperature range $`1.5T/T_c3`$, the spectral density is peaked at small energies. This peak was next analyzed using a parametrised form of the Bayesian prior and the electrical conductivity of the plasma was extracted. This was then used to predict the soft photon emissivity of the QCD plasma.

Rough estimates of typical transport related time and length scales in the QCD plasma can be obtained by using the extracted value of the electrical conductivity in conjunction with the simple transport formula $`\sigma =C_{EM}n_q\tau _q/m`$ ($`C_{EM}`$ is nothing but the average charge square: $`e^2`$). If the number density of quarks, $`n_q`$, is substituted by the corresponding entropy density, and the screening mass used for $`m`$, then one finds the quark mean free time $`\tau _q0.3`$ fm. This is also the time scale for the persistence of charge, isospin, strangeness and baryon number fluctuations in the plasma amy2 . A possible experimental check could be to determine the mean free path of very soft off-shell photons ($`\omega 0.15(2\pi T)200`$ MeV). These are only 20 times longer than $`\tau _q`$, i.e., about 6 fm. The fireball at RHIC may be marginally transparent to such photons, but at LHC the fireball size could be large enough to attenuate the intensity of very soft photons. Such an observation would constitute direct evidence for short mean free paths in the plasma.

Some estimates of other transport coefficients can be obtained if one assumes that the mean free time of gluons is $`\tau _q/2`$, since they should be related by colour factors. Then simple transport formulæ treated in the same approximation as before lead us to the estimate that the dimensionless ratio $`\eta /S0.2`$, where $`S`$ is the entropy density of the plasma and $`\eta `$ is the shear viscosity. This ratio is of the same order of magnitude as extracted from present heavy-ion data teany  and obeys a bound conjectured in son . It would be useful to make a direct measurement of the shear viscosity on the lattice.

Many interesting lines of research are relegated to the future. The dilepton emissivity away from $`\omega =0`$ is a conceptually simple extension, but requires further numerical work, as explained before. Extending these measurements closer toward $`T_c`$ where correlation lengths grow larger saumen  is of obvious importance, but outside the scope of this paper, as is the extension to dynamical QCD. The interesting question of the effectiveness of linear response theory, and hence of the Kubo formulae closer to $`T_c`$, can perhaps be probed using the non-linear susceptibilities defined in zakopane .


text:
# Wavelet analysis of epileptic spikes

(January 17, 2003)

## Abstract

Interictal spikes and sharp waves in human EEG are characteristic signatures of epilepsy. These potentials originate as a result of synchronous, pathological discharge of many neurons. The reliable detection of such potentials has been the long standing problem in EEG analysis, especially after long-term monitoring became common in investigation of epileptic patients. The traditional definition of a spike is based on its amplitude, duration, sharpness, and emergence from its background. However, spike detection systems built solely around this definition are not reliable due to the presence of numerous transients and artifacts. We use wavelet transform to analyze the properties of EEG manifestations of epilepsy. We demonstrate that the behavior of wavelet transform of epileptic spikes across scales can constitute the foundation of a relatively simple yet effective detection algorithm.

Recordings of human brain electrical activity (EEG) have been the fundamental tool for studying the dynamics of cortical neurons since 1929. Even though the gist of this technique has essentially remained the same, the methods of EEG data analysis have profoundly evolved during the last two decades. In 1985 Babloyantz et al. demonstrated that certain nonlinear measures, first introduced in the context of chaotic dynamical systems, changed during slow-wave sleep Babloyantz et al. (1985). The flurry of research work that followed this discovery focused on the application of nonlinear dynamics in quantifying brain electrical activity during different mental states, sleep stages, and under the influence of the epileptic process (for a review see for example West et al. (1995); Pritchard and Duke (1992)). It must be emphasized that a straightforward interpretation of neural dynamics in terms of such nonlinear measures as the largest Lyapunov exponent or the correlation dimension is not possible since most biological time series, such as EEG, are nonstationary and consequently do not satisfy the assumptions of the underlying theory. On the other hand, traditional power spectral methods are also based on quite restrictive assumptions but nevertheless have turned out to be successful in some areas of EEG analysis. Despite these technical difficulties, the number of applications of nonlinear time series analysis has been growing steadily and now includes the characterization of encephalopaties Stam et al. (1999), monitoring of anesthesia depth Widman et al. (2000), characteristics of seizure activity Casdagli et al. (1997), and prediction of epileptic seizures Lehnertz and Elger (1998). Several other approaches are also used to elucidate the nature of electrical activity of the human brain ranging from coherence measures Nunez et al. (1997, 1999) and methods of nonequilibrium statistical mechanics Ingber (1982) to complexity measures Rezek and Roberts (1998); van Putten and Stam (2001).

One of the most important challenges of EEG analysis has been the quantification of the manifestations of epilepsy. The main goal is to establish a correlation between the EEG and clinical or pharmacological conditions. One of the possible approaches is based on the properties of the interictal EEG (electrical activity measured between seizures) which typically consists of linear stochastic background fluctuations interspersed with transient nonlinear spikes or sharp waves. These transient potentials originate as a result of a simultaneous pathological discharge of neurons within a volume of at least several mm3.

The traditional definition of a spike is based on its amplitude, duration, sharpness and emergence from its background Gotman and Gloor (1976); Gotman and Wang (1982). However, automatic epileptic spike detection systems based on this direct approach suffer from false detections in the presence of numerous types of artifacts and non-epileptic transients. This shortcoming is particularly acute for long-term EEG monitoring of epileptic patients which became common in 1980s. To reduce false detections Gotman and Wang Gotman and Wang (1991) made the process of spike identification dependent upon the state of EEG (active wakefulness, quiet wakefulness, desynchronized EEG, phasic EEG, and slow EEG). This modification leads to significant overall improvement provided that state classification is correct.

Diambra and Malta Diambra and Malta (1999) adopted nonlinear prediction for epileptic spike detection. They demonstrated that when the model’s parameters are adjusted during the “learning”phase to assure good predictive performance for stochastic background fluctuations, the appearance of an interictal spike is marked by a very large forecasting error. This novel approach is appealing because it makes use of changes in EEG dynamics. One expects good nonlinear predictive performance when the dynamics of the EEG interval used for building up the model is similar to the dynamics of the interval used for testing. However, it is uncertain at this point whether it is possible to develop a robust spike detection algorithm based solely on this idea.

As Clark et al. put it succinctly, automatic EEG analysis is a formidable task because of the lack of “…features that reflect the relevant information ”Clark et al. (1995). Another difficulty is the nonstationary nature of the spikes and the background in which they are embedded. One technique developed for the treatment of such nonstationary time series is wavelet analysis Mallat (1998); Unser and Aldroubi (1996). The goal of this paper is to characterize the epileptic spikes and sharp waves in terms of the properties of their wavelet transforms. In particular, we search for features which could be important in the detection of epileptic events.

The wavelet transform is an integral transform for which the set of basis functions, known as wavelets, are well localized both in time and frequency. Moreover, the wavelet basis can be constructed from a single function $`\psi (t)`$ by means of translation and dilation:

$$\psi _{a;t_0}=\psi \left(\frac{tt_0}{a}\right).$$
(1)
$`\psi (t)`$ is commonly referred to as the mother function or analyzing wavelet. The wavelet transform of function $`h(t)`$ is defined as

$$W(a,t_0)=\frac{1}{\sqrt{a}}_{\mathrm{}}^{\mathrm{}}h(t)\psi _{a;t_0}^{}𝑑t,$$
(2)
where $`\psi ^{}(t)`$ denotes the complex conjugate of $`\psi (t)`$. The continuous wavelet transform of a discrete time series $`\{h_i\}_{i=0}^{N1}`$ of length $`N`$ and equal spacing $`\delta t`$ is defined as

$$W_n(a)=\sqrt{\frac{\delta t}{a}}\underset{n^{}=0}{\overset{N1}{}}h_n^{}\psi ^{}\left[\frac{(n^{}n)\delta t}{a}\right].$$
(3)
The above convolution can be evaluated for any of $`N`$ values of the time index $`n`$. However, by choosing all $`N`$ successive time index values, the convolution theorem allows us to calculate all $`N`$ convolutions simultaneously in Fourier space using a discrete Fourier transform (DFT). The DFT of $`\{h_i\}_{i=0}^{N1}`$ is

$$\widehat{h}_k=\frac{1}{N}\underset{n=0}{\overset{N1}{}}h_ne^{2\pi ikn/N},$$
(4)
where $`k=0,\mathrm{},N1`$ is the frequency index. If one notes that the Fourier transform of a function $`\psi (t/a)`$ is $`|a|\widehat{\psi }(af)`$ then by the convolution theorem

$$W_n(a)=\sqrt{a\delta t}\underset{k=0}{\overset{N1}{}}\widehat{h}_n\psi ^{}(af_k)e^{2\pi if_kn\delta t},$$
(5)
frequencies $`f_k`$ are defined in the conventional way. Using (5) and a standard fast Fourier transform (FFT) routine it is possible to efficiently calculate the continuous wavelet transform (for a given scale $`a`$) at all $`n`$ simultaneously Torrence and Compo (1998). It should be emphasized that formally equation (5) does not yield the discrete linear convolution corresponding to (3) but rather a discrete circular convolution in which the shift $`n^{}n`$ is taken modulo $`N`$. However, in the context of this work, this problem does not give rise to any numerical difficulties. This is because, for purely practical reasons, the beginning and the end of the analyzed part of data stream are not taken into account during the EEG spike detection.

From a plethora of available mother wavelets, we employ the Mexican hat

$$\psi (t)=\frac{2}{\sqrt{3}}\pi ^{1/4}(1t^2)e^{t^2/2}$$
(6)
which is particularly suitable for studying epileptic events.

In the top panel of Fig. 1 we present two pieces of the EEG recording joined at approximately $`t=1s`$. The digital 19 channel recording sampled at 240 Hz was obtained from a juvenile epileptic patient according to the international 10-20 standard with the reference average electrode. The epileptic spike in this figure (marked by the arrow) is followed by two artifacts. The bottom panel of Fig. 1 displays the contour map of the absolute value of Mexican hat wavelet coefficients $`W(a,t_0)`$. It is apparent that the red prominent ridges correspond to the position of either spike or the motion artifacts. What is most important, for small scales, $`a`$, the values of the wavelet coefficients for the spike’s ridge are much larger than those for the artifacts. The peak value along the spike ridge corresponds to $`a=7`$. In sharp contrast, for the range of scales used in Fig. 1 the absolute value of coefficients $`W(a,t_0)`$ for the artifacts grow monotonically with $`a`$.

The question arises as to whether the behavior of the wavelet transform as a function of scale can be used to develop a reliable detection algorithm. The first step in this direction is to use the normalized wavelet power

$$w(a,t_0)=W^2(a,t_0)/\sigma ^2$$
(7)
instead of the wavelet coefficients to reduce the dependence on the amplitude of the EEG recording. In the above formula $`\sigma ^2`$ is the variance of the portion of the signal being analyzed (typically we use pieces of length 1024 for EEG tracings sampled at 240 Hz). In actual numerical calculations we prefer to use the square of $`w(a,t_0)`$ to merely increase the range of values analyzed during the spike detection process. In Fig. 2 $`w^2`$ for the signal used in Fig. 1 is plotted for three scales $`A=3`$, $`B=7`$ and $`C=20`$.

In the most straightforward approach, we identify an EEG transient potential as a simple or isolated epileptic spike if and only if:

$``$ the value of $`w^2`$ at $`a=7`$ is greater than a predetermined threshold value $`T_1`$,

$``$ the square of normalized wavelet power decreases from scale $`a=7`$ to $`a=20`$,

$``$ the value of $`w^2`$ at $`a=3`$ is greater than a predetermined threshold value $`T_2`$.

The threshold values $`T_1`$ and $`T_2`$ may be considered as the model’s parameters which can be adjusted to achieve the desired sensitivity (the ratio of detected epileptic events to the total number of epileptic events present in the analyzed EEG tracing) and selectivity (the ratio of epileptic events to the total number of events marked by the algorithm as epileptic spikes).

While this simple algorithm is quite effective for simple spikes such as one shown in Fig. 1 it fails for the common case of an epileptic spike accompanied by a slow wave with comparable amplitude. The example of such complex is given in Fig. 3(a). The overlap of the negative tail of the Mexican hat with the slow wave yields the inherently low values of $`w^2`$ at scale $`A`$ (panel (b)) and scale $`B`$ (panel (c)) as compared to those characteristic of the “isolated”spike. Nevertheless, the normalized wavelet power does decrease from scale $`B`$ to $`C`$. Consequently, in the same vein as the argument we presented above, we can develop an algorithm which detects the epileptic spike in the vicinity of a slow wave by calculating the following linear combination of wavelet transforms:

$$\stackrel{~}{W}(a,t_0)=c_1W(a,t_0)+c_2W(a_s,t_0+\tau )$$
(8)
and checking whether the square of corresponding normalized power $`\stackrel{~}{w}(a,t_0)=\stackrel{~}{W}^2(a,t_0)/\sigma ^2`$at scales $`a=3`$ and $`a=7`$ exceeds the threshold value $`\stackrel{~}{T}_1`$ and $`\stackrel{~}{T}_2`$, respectively. The second term in (8) allows us to detect the slow wave which follows the spike. The parameters $`a_s`$ and $`\tau `$ are chosen to maximize the overlap of the wavelet with the slow wave. For the Mexican hat we use $`a_s=28`$ and $`\tau =0.125s`$. By varying the values of coefficients $`c_1`$ and $`c_2`$, it is possible to control the relative contribution of the spike and the slow wave to the linear combination (8).

For testing purposes, we built up the database of artifacts and spikes. We made available some of these EEG tracings at hom  along with the examples of the numerical calculations. While the analysis of the pieces of EEG recordings such as those shown in Fig. 2 and 3 is essential in determining the generic properties of epileptic events, it can hardly reflect the difficulties one can encounter in interpretation of clinical EEG. Therefore we selected four challenging EEG tracings with 340 epileptic events. The algorithm described in this work marked 356 events out of which 239 turned out to be the epileptic events. Thus the sensitivity of the algorithm was 70% and its selectivity was equal to 67%. We then analyzed the same tracings with the leading commercial spike detector developed by the Persyst Development Corporation (Insight 2001.07.12). This software marked 654 events out of which 268 were epileptic events. Thus slightly better sensitivity of 79% was achieved at the expense of the low 41% selectivity. The performance of preliminary numerical implementation of the detection algorithm presented in this work is excellent and allows to process 24 hour EEG recording (19 channels sampled at 240 Hz) in a matter of minutes on the average personal computer.

The goal of wavelet analysis of the two types of spikes, presented in this paper, was to elucidate the approach to epileptic events detection which explicitly hinges on the behavior of wavelet power spectrum of EEG signal across scales and not merely on its values. Thus, this approach is distinct not only from the detection algorithms based upon discrete multiresolution representations of EEG recordings Blanco et al. (1996); D’Attellis et al. (1997); Sartoretto and Ermani (1999); Calvagno et al. (2000); Gutiérrez et al. (2001) but also from the method developed by Senhadji and Wendling which employs continuous wavelet transform Senhadji and Wendling (2002).

Epilepsy is a common disease which affects 1-2% of the population and about 4% of children Jallon (2002). In some epilepsy syndromes interictal paroxysmal discharges of cerebral neurons reflect the severity of the epileptic disorder and themselves are believed to contribute to the progressive disturbances in cerebral functions (eg. speech impairment, behavioral disturbances) Engel (2001). In such cases precise quantitative spike analysis would be extremely important. The epileptic event detector described in this paper was developed with this particular goal in mind and its application to the studies of the Landau-Kleffner syndrome will be presented elsewhere.


text:
# The QCD Pomeron in Ultraperipheral Heavy Ion Collisions: II. Heavy Quark Production

## I INTRODUCTION

The behavior of $`ep/pp`$ scattering in the limit of high center-of-mass energy $`\sqrt{s}`$ and fixed momentum transfer is one of the outstanding open questions in the theory of the strong interactions. In the late 1970s, Lipatov and collaborators  established the papers which form the core of our knowledge of Regge limit (high energy limit) of Quantum Chromodynamics (QCD). The physical effect that they describe is often referred to as the QCD Pomeron, with the evolution described by the BFKL equation. The simplest process where this equation applies is the high energy scattering between two heavy quark-antiquark states, i.e. the onium-onium scattering. For a sufficiently heavy onium state, high energy scattering is a perturbative process since the onium radius gives the essential scale at which the running coupling $`\alpha _s`$ is evaluated. This process was proposed as a gedanken experiment to investigate the high energy regime of QCD in Refs.  (see also Refs. ). In the dipole picture , the heavy quark-antiquark pair and the soft gluons in the limit of large number of colors $`N_c`$ are viewed as a collection of color dipoles. In this case, the cross section can be understood as a product of the number of dipoles in one onium state, the number of dipoles in the other onium state and the basic cross section for dipole-dipole scattering due to two-gluon exchange. At leading order (LO), the cross section grows rapidly with the energy ($`\sigma \alpha _s^2e^{(\alpha _{IP}1)Y}`$, where $`(\alpha _{IP}1)=\frac{4\alpha _sN_c}{\pi }\mathrm{ln}2`$ and $`Y=\mathrm{ln}s/Q^2`$) because the number of dipoles in the light cone wave function grows rapidly with the energy. Therefore, in principle, the analysis of the energy dependence of the cross section could disentangle the QCD Pomeron effects. However, for $`ep/pp`$ colliders, the study of the QCD Pomeron is made difficult by the fact that the cross section is influenced by both short and long distance physics. Only when specific conditions are satisfied is that one can expect to determine the QCD Pomeron effects. Some examples are the forward jet production in deeply inelastic events at low values of the Bjorken variable $`x`$ in lepton-hadron scattering and jet production at large rapidity separations in hadron-hadron collisions.

Over the past few years much experimental effort has been devoted towards to observe the QCD Pomeron in those processes, but the conclusions are not unambiguous (see e.g. Ref. ). This fact is mainly associated with the theoretical uncertainty present in the high energy limit of QCD. Only recently, after an effort of ten years, the next-to-leading order (NLO) corrections to the BFKL equation were obtained (for a review on NLO BFKL corrections, see e.g. Ref.  and references therein). Moreover, since the BFKL equation predicts that for $`s\mathrm{}`$ the corresponding cross section rises with a power law of the energy, violating the Froissart bound, new dynamical effects associated with the unitarity corrections are expected to stop further growth of the cross sections . This expectation can be easily understood: while for large transverse momentum $`k_{}`$, the BFKL equation predicts that the mechanism $`ggg`$ populates the transverse space with a large number of small size gluons per unit of rapidity (the transverse size a gluon with momentum $`k_{}`$ is proportional to $`1/k_{}^2`$), for small $`k_{}`$ the produced gluons overlap and fusion processes, $`ggg`$, are equally important. Considering this process, the rise of the gluon distribution with transverse momenta below a typical scale, energy dependent, called saturation scale $`Q_s`$ is reduced, restoring the unitarity. It is important to salient that Golec-Biernat and Wusthoff  have shown that a saturation model is able to describe the DESY $`ep`$ collider HERA data, in particular the transition from the perturbative to the nonperturbative photoproduction region, including the inclusive and diffractive cross sections. Therefore, the description of the QCD Pomeron still is an open question, which deserves more detailed analysis.

This situation should be improved in the future with the next generation of linear colliders (For a recent review of BFKL searches, see e.g. Ref. ). In particular, future $`e^+e^{}`$ colliders probably will allow to discriminate between BFKL and saturation predictions  (For other analysis of photon-photon collisions see Refs. ). Such reaction presents analogies with the process of scattering of two onia discussed above. Although the onium-onium scattering is a gedanken experiment, off-shell photon scattering at high energy in $`e^+e^{}`$ colliders, where the photons are produced from the leptons beams by bremsstrahlung, plays a similar role. In these two-photon reactions, the photon virtualities can be made large enough to ensure the applicability of the perturbative methods. Similarly, we can expect to test the QCD Pomeron in heavy flavor production in two-photon collisions, where the hard scale is provided by the heavy quark mass.

In the last years many authors have studied in detail the referred process. From the view of the color dipole picture, in Ref.  the contribution of the dipole-dipole interaction was taken into account considering a generalization of the original saturation model  to the two-photon reactions. There, the results present good agreement with data from $`e^+e^{}`$ colliders, remaining some room to discuss issues concerning resolved contributions. In Ref. , these questions were revisited using the saturation model, where threshold effects were addressed in an accurate way as well as the hadronic single-resolved contribution was computed within the dipole formalism. A current open question is the discrepancy between the theoretical results and the data on heavy quark production, mostly for the bottom case, in the $`\gamma p`$ photoproduction and two-photon reactions. This feature is also corroborated by the dipole-dipole approach from Ref. . As found in Ref. , even the inclusion of resolved contributions does not remove completely the observed deficit. Quite promising results are obtained in Ref. , where heavy quark production in two photon reaction is calculated within the $`k_{}`$-factorization formalism and the unintegrated gluon density for the photon is determined using the saturation model.

¿From the view of the BFKL approach, there are several calculations using the leading logarithmic approximation  and considering the next-to-leading corrections to the total cross section $`\gamma ^{}\gamma ^{}`$ process . In particular, a stable next-to-leading order program relying on the BLM optimal scale setting  has produced good results with OPAL and L3 data at LEP2 . However, the heavy quark production from real two-photon process is still an open question. Only the effects of charm mass in the $`\gamma ^{}\gamma ^{}`$ total cross section for the kinematical region of L3 at LEP was analyzed by Bartels et al.  considering the LO BFKL dynamics. We address that issue in our analysis, considering a BFKL dipole approximation to the QCD dynamics. Similarly to the previous analysis on two-photon physics, one shows that the next generation of colliders will be able to discriminate the dynamics associated to the QCD Pomeron.

Recently, we have proposed in Ref.  to investigate QCD Pomeron effects in a different context, namely in photon-photon scattering at ultraperipheral heavy ion collisions. In this case, the cross sections are enhanced since the $`\gamma \gamma `$ luminosity increases as $`Z^4`$, where $`Z`$ is the atomic number . There, we have analyzed the double diffractive $`J/\mathrm{\Psi }`$ production in $`\gamma \gamma `$ collisions, with the photons coming from the Weizsäcker - Williams spectrum of the nuclei. For that process our results have indicated that future experimental analysis can be useful to discriminate the QCD dynamics at high energies. Here we extend our proposal for heavy quark production. In principle, at ultraperipheral heavy ion collisions, a photon stemming from the electromagnetic field of one of the two colliding nuclei can interact with one photon of the other nucleus (two-photon process) or can penetrate into the other nucleus and interact with its hadrons (photon-nucleus process). Therefore, heavy quarks can be produced in photo-nucleus interactions as well as photon-photon interactions. The first case has been extensively discussed in literature  (For a review, see Ref. ) and recently revised in different approaches . On the other hand, heavy quark production in two-photon collisions has been discussed essentially in Refs. . In these analysis, the underlying contribution comes from the QED quark-box diagram (QPM) at LO \[Fig. 1 (top panel)\], which dominates in the available energies at RHIC, $`W_{\gamma \gamma }6`$ GeV. In Ref. , besides the direct contribution also the resolved contributions were carefully calculated and the latter is shown to be negligible. However, at the LHC experiment, where the two-photon energies can reach to $`W_{\gamma \gamma }200`$ GeV, the QCD diagrams become important and the QCD Pomeron should dominate in this limit \[Fig. 1 (bottom panel)\]. Consequently, an enhancement of the cross section associated to the dynamics is expected. We will address these issues in the section III, complementing the detailed study presented in Ref. .

Here, we will restrict our analysis for the two-photon process and its potentiality to investigate the QCD dynamics. Similarly as in Ref. , our goal here is twofold: to analyze the potentiality of this process to constraint the QCD dynamics at high values of energy and to provide reliable estimates for the cross sections concerning that reaction. Relativistic heavy-ion collisions are a potentially prolific source of $`\gamma \gamma `$ collisions at high energy colliders. The advantage of using heavy ions is that the cross sections varies as $`Z^4\alpha ^4`$ rather then only $`\alpha ^4`$. Moreover, the maximum $`\gamma \gamma `$ collision energy $`W_{\gamma \gamma }`$ is $`2\gamma /R_A`$, about 6 GeV at RHIC and 200 GeV at LHC, where $`R_A`$ is the nuclear radius and $`\gamma `$ is the center-of-mass system Lorentz factor of each ion. For two-photon collisions, the cross section for the reaction $`AAAAQ\overline{Q}`$, where $`Q`$ $`=c,b`$, will be given by

$`\sigma _{AAAAQ\overline{Q}}(s)={\displaystyle 𝑑\tau \frac{d_{\gamma \gamma }}{d\tau }\widehat{\sigma }_{\gamma \gamma Q\overline{Q}}(\widehat{s})},`$ (1)
where $`\tau =\widehat{s}/s`$, $`\widehat{s}=W_{\gamma \gamma }^2`$ is the square of the center of mass (c.m.s.) system energy of the two photons, $`s`$ of the ion-ion system, $`d_{\gamma \gamma }/d\tau `$ is the photon luminosity and $`\widehat{\sigma }_{\gamma \gamma Q\overline{Q}}(\widehat{s})`$ is the cross section of the $`\gamma \gamma `$ interaction (For details related to the numerical expressions, see Ref. ). Our approach excludes possible final state interactions of the produced particles with the colliding particles, allowing reliable calculations of ultraperipheral heavy ion collisions. Therefore, to estimate the heavy quark production it is only necessary to consider a suitable QCD model for photon-photon interaction.

This paper is organized as follows. In Sec. II, we calculate the total cross sections for the process $`\gamma \gamma Q\overline{Q}`$ based on the color dipole picture, using the dipole cross section given by the saturation model and BFKL approach. In Sec. III, the corresponding heavy quark cross section in coherent ultraperipheral heavy ion collisions are presented, with a particular emphasis in the LHC energies. In the last section we drawn our conclusions and summarize the main results.

## II Heavy-quark production at the photon level

The heavy quark production via photon-photon collision is a powerful process testing the applicability of perturbative QCD methods. The photon induced reactions have the advantage of being clear concerning the probing structure and the large mass of the heavy quarks produced allows perturbative physics to be employed. The mechanisms of production of such process have been systematically studied, despite the available data being modest and having large uncertainties (see Ref.  and references therein). An important experimental verification is that the data on open heavy quark pair production, mostly for bottom quark, stay above theoretical calculations by a sizeable factor. Such a feature has generated several studies which rely on higher orders in the collinear factorization approach or on the successful results based on the semihard approach.

An important approach describing two-photon process is given by the color dipole formalism . The simple physical picture is provided by the photon splitting in a quark-antiquark pair (virtual components in the transverse plane), called color dipole, long after the interaction and the further scattering of these pairs with the target. Such an approach has produced a consistent and unified description of high energy photon induced processes in deep inelastic scattering and two-photon reactions, including the deep inelastic diffractive dissociation. In particular, in the two-photon interactions, the target is also given by a color dipole. The wave functions describing the photon virtual fluctuation are completely determined from perturbative methods, whereas the effective color dipole cross section has to be modelled. It may contain nonperturbative and higher twist contributions as well as parton saturation effects as we should discuss latter on.

The $`Q\overline{Q}`$ pair production mechanism can be separated into three main contributions: (i) the quark box diagram (QPM), with photons coupling to the same heavy quark line, which dominates at not so large photon-photon c.m.s. energies, $`W_{\gamma \gamma }`$ \[Fig. 1 (top panel)\]; (ii) the dipole-dipole interaction, where both photons fluctuate in quark-antiquark pairs and interact through gluonic exchange, which dominates at high energies \[Fig. 1 (bottom panel)\]; (iii) single and double resolved contributions, where one of the photons fluctuates into the vector mesons ($`\rho `$, $`\omega `$, $`\varphi `$, $`\mathrm{}`$) interacting with the remaining dipole, being a background to the dipole-dipole reaction. The first contribution is well known, providing a $`W_{\gamma \gamma }^2`$ behavior modulo logarithmic corrections, given by for example on Refs. . The second one will be reproduced in the following, making use of a sound model for the dipole cross section . The third one will be disregarded, but we quote the study presented in Ref. , where these pieces were calculated for the first time in the color dipole approach.

The color dipole formulation gives the following expression for the heavy quark pair production in a (real) two-photon reactions at high energies,

$`\sigma ^{\mathrm{dd}}(\gamma \gamma Q\overline{Q})`$ $`=`$ $`{\displaystyle \underset{q_1Q}{}}{\displaystyle |\mathrm{\Psi }_{q_1\overline{q}_1}(𝒓_1,z_1)|^2|\mathrm{\Psi }_{Q\overline{Q}}(𝒓_2,z_2)|^2\sigma _{\mathrm{dd}}(𝒓_1,𝒓_2,\stackrel{~}{x}_{ab})d^2𝒓_1d^2𝒓_2𝑑z_1𝑑z_2}`$ (2)
$`+`$ $`{\displaystyle \underset{q_2Q}{}}{\displaystyle |\mathrm{\Psi }_{Q\overline{Q}}(𝒓_1,z_1)|^2|\mathrm{\Psi }_{q_2\overline{q}_2}(𝒓_2,z_2)|^2\sigma _{\mathrm{dd}}(𝒓_1,𝒓_2,\stackrel{~}{x}_{ab})d^2𝒓_1d^2𝒓_2𝑑z_1𝑑z_2},`$ (3)
where $`\mathrm{\Psi }_{q\overline{q},(Q\overline{Q})}(𝒓,z)`$ are the light (heavy) quark-antiquark wave functions of the photon in the mixed representation. The transverse separation (dipole size) of the respective quark pair is denoted by $`𝒓`$ and their longitudinal momentum fraction denoted by $`z`$. The dipole-dipole cross section, $`\sigma ^{\mathrm{dd}}`$, depends on the dipole sizes and on the parameter $`\stackrel{~}{x}_{ab}`$, driving the energy behavior.

First, we write down the expressions for the wave functions squared including the respective quark helicities $`T,L`$ and photon virtuality defined as $`Q^2`$,

$`|\mathrm{\Psi }_T(𝒓,z)|^2`$ $`=`$ $`{\displaystyle \frac{6\alpha _{em}}{4\pi ^2}}{\displaystyle \underset{f}{}}e_f^2[z^2+(1z)^2]\epsilon _f^2K_1^2(\epsilon _fr)+m_f^2K_0^2(\epsilon _fr),`$ (4)
$`|\mathrm{\Psi }_L(𝒓,z)|^2`$ $`=`$ $`{\displaystyle \frac{6\alpha _{em}}{4\pi ^2}}{\displaystyle \underset{f}{}}e_f^2\mathrm{\hspace{0.17em}4}Q^2z^2(1z)^2K_0^2(\epsilon _fr),`$ (5)
where $`e_f`$ and $`m_f`$ stands for the charge and mass of the quark having flavor $`f`$ and $`K_{0,1}`$ are the McDonald-Bessel functions. The auxiliary variable $`\epsilon `$ is defined as $`\epsilon ^2=z(1z)Q^2+m_f^2`$ which in the real photon case, $`Q^2=0`$, considered here, simplifies to $`\epsilon =m_f`$. Also in such case, the longitudinal photon component is obviously suppressed. The light quark masses are considered to be equal and an effective value $`m_q0.2`$, is determined from fitting two-photon data .

Lets discuss now a specific model for the dipole-dipole cross section. Here, we follows the formulation of the saturation model applied to two-photon reactions . The cross section describing the dipole interaction satisfies the saturation property, i.e. a constant dependence on energy at large dipole sizes, and color transparency behavior at small dipole sizes. The interface between this two regions is controlled by the saturation scale, which is energy dependent and has been determined from fitting data on small-$`x`$ deep inelastic data . Explicitly, the dipole-dipole cross section takes the eikonal-like form,

$`\sigma _{\mathrm{dd}}(𝒓_1,𝒓_2,\stackrel{~}{x}_{ab})=\stackrel{~}{\sigma }_0\left[1\mathrm{exp}\left({\displaystyle \frac{\overline{r}^2}{4R_0^2(\stackrel{~}{x}_{ab})}}\right)\right],`$ (6)
with the overall normalization given by $`\stackrel{~}{\sigma }_0=(2/3)\sigma _0`$ inspired in the quark counting rule and where $`\sigma _0`$ is determined for the proton case analysis in Ref. . The effective radius $`\overline{r}`$ is constructed is such way to reproduce the dipole-hadron parameterization, i.e. $`\overline{𝒓}^2𝒓_1^2(𝒓_2^2)`$ in the dipole configurations $`𝒓_2^2𝒓_1^2(𝒓_1^2𝒓_2^2)`$. Here, we will consider two possibilities giving good results as found in Ref. : (I) $`\overline{r}^2=r_1^2r_2^2/(r_1^2+r_2^2)`$, preferable against two-photon data; (II) $`\overline{r}^2=\mathrm{min}(r_1^2,r_2^2)`$, preferable for charm quark pair production data. The fitted mass for the light quarks takes the value $`m_q=0.21`$ GeV considering the effective radius (I) and $`m_q=0.23`$ GeV using radius (II). The saturation radius $`R_0`$ squared reads as,

$`R_0^2(\stackrel{~}{x}_{ab})=\left({\displaystyle \frac{\stackrel{~}{x}_{ab}}{x_0}}\right)^\lambda \mathrm{GeV}^2,`$ (7)
where the parameters $`x_0`$ and $`\lambda `$ are also determined in the original saturation model in the proton case . The saturation scale, $`Q_s(\stackrel{~}{x}_{ab})`$, is obtained by the inverse of the saturation radius described above. The effective variable $`\stackrel{~}{x}`$ depends on the two-photon energy and the quark masses and reads as,

$`\stackrel{~}{x}_{ab}={\displaystyle \frac{4m_a^2+4m_b^2}{W^2}},`$ (8)
where in the heavy quark production $`a`$ and $`b`$ corresponds to a light and a heavy (or vice versa) quark. For low energy (large $`\stackrel{~}{x}`$), the dipole cross section needs to be supplemented by threshold effects. Following , they are taken into account through the multiplicative factor $`(1\stackrel{~}{x}_{ab})^5`$, coming from dimensional-counting rules. However, a different procedure is adopted in Ref. , where the threshold is considered by imposing the kinematical constraint $`M_{q\overline{q}}+M_{Q\overline{Q}}<W_{\gamma \gamma }`$ in the integrations at Eq. (3), where $`M`$ is the invariant mass of the pair. There, also the variable $`\stackrel{~}{x}`$ takes a different form than Eq. (8), depending explicitly on $`z_1`$ and $`z_2`$.

In Fig. 2 we reproduce the results on charm production considering the saturation model (SAT-MOD), regarding the quark box contribution (QPM) and the dipole-dipole approach using the saturation model. The theoretical predictions are compared with the experimental results from the L3 Collaboration . It should be stressed that results for the bottom case are straightforward, changing only the quark charge and mass, $`e_b=1/3`$ and $`m_b=4.5`$ GeV. The charm mass is taken as $`m_c=1.3`$ GeV and the two prescriptions for the effective radius are presented. The QPM contribution (long-dashed curve) is very important at low energy, starting to be negligible around $`W_{\gamma \gamma }=10`$ GeV. Therefore, this contribution should dominate the ultraperipheral heavy quark production at RHIC energies, where the maximal $`W_{\gamma \gamma }`$ reached is of order of few GeV’s. However, above $`W_{\gamma \gamma }10`$ GeV, the dipole-dipole approach dominates the physical description, providing an increase of the total cross section. The results using the effective radius (I) (solid line) and (II) (dashed line) are presented. The contribution from dipole-dipole interaction is shown separately for the radius (I), represented by the dot-dashed line. This contribution plays a very important role at LHC energies, where the two-photon energy in ultraperipheral collisions can reach up to $`W160`$ GeV.

Regarding the resolved (single and double) contributions, we have not considered them here. We quote Ref. , where these contributions are calculated using the saturation model within the dipole-dipole picture. There, it was found that for the charm case the single resolved piece constitutes about 30 % of the main contribution. In the bottom case, these contributions are quite smaller.

A prominent challenge in high energy physics concerns the issue what is the correct dynamics driving the interaction at asymptotic energies. It is basically dominated by gluonic interaction and the saturation model discussed above successfully interpolates between the soft (non-perturbative) and hard (perturbative) regimes, whose transition relies on the parton saturation phenomenon. However, the physics of two-photon collisions is an optimal tool testing the BFKL approach . In its dipole formulation, the BFKL dynamics is suitable to describe the interaction of small size color dipoles from the interacting photons. A lot of work have been done using the LO BFKL approach and also attempts to include NLO corrections in two-photon reactions . However, heavy quark production in real photon-photon collisions considering the BFKL approach was not addressed in the current literature. Below, we present our results for this case.

To estimate the heavy quark production in the BFKL approach we need a model for the dipole-dipole cross section. In particular, the BFKL dipole cross section can be obtained from the momentum representation and it can be considered as a solution of the linearized Balitsky-Kovchegov equation in the coordinate space . The LO BFKL solution is completely known in the momentum space, providing a steep energy dependence associated with the hard Pomeron intercept $`1+\omega _{IP}`$ and the typical diffusion feature. For our purpose here, we will consider the following expression for the dipole-dipole cross section from a BFKL dynamics,

$`\sigma _{\mathrm{dd}}(𝒓_1,𝒓_2,\stackrel{~}{x}_{ab})={\displaystyle \frac{4\pi ^2}{3}}𝒞\left({\displaystyle \frac{\overline{r}^2}{\overline{r}_0^2}}\right)^{\frac{1}{2}}{\displaystyle \frac{e^{\omega \overline{\alpha }_s\mathrm{ln}(1/\stackrel{~}{x}_{ab})}}{\sqrt{2\pi \stackrel{~}{\omega }\overline{\alpha }_s\mathrm{ln}(1/\stackrel{~}{x}_{ab})}}}\mathrm{exp}\left[{\displaystyle \frac{\mathrm{ln}^2(\overline{r}^2/\overline{r}_0^2)}{2\stackrel{~}{\omega }\overline{\alpha }_s\mathrm{ln}(1/\stackrel{~}{x}_{ab})}}\right](1\stackrel{~}{x}_{ab})^5,`$ (9)
where $`\overline{\alpha }_s=N_c\alpha _s/\pi `$, with number of color $`N_c=3`$. Moreover, for light quarks, $`𝒞=1.0`$, $`\overline{r}_0^2=1/\mathrm{\Lambda }_{\mathrm{QCD}}^2`$, where one has used $`\mathrm{\Lambda }_{QCD}=0.23`$. Still, $`\omega 2\psi (1)2\psi (1/2)=4\mathrm{ln}\mathrm{\hspace{0.17em}2}`$ and $`\stackrel{~}{\omega }=28\xi (3)`$. The Pomeron intercept is given by $`\omega _{IP}=1+\overline{\alpha }_s\omega `$. The prescription used here is not unique. For instance, for dipole-dipole collisions, the choice would be $`\overline{r}=\mathrm{min}(𝒓_1,𝒓_2)`$ and $`\overline{r}_0=\mathrm{max}(𝒓_1,𝒓_2)`$, having impact on the final results.

In Fig. 3 are presented the results using the dipole-dipole cross section of Eq. (9), considering the effective radius (I) and taking into account three different values for the Pomeron intercept. The QPM contribution was also included. In our further analysis one considers the intermediate value $`\omega _{IP}=1.47`$, which is steeper than the results from the saturation model. It should be stressed that further phenomenology can be done, for instance taking the overall normalization, the scale $`\overline{r}_0`$ and the Pomeron intercept $`\omega _{IP}`$ as free parameters. The differences between the predictions of the saturation model and the BFKL approach are already sizeable in the energy interval presented here, with important implications in the future TESLA experiment.

As a final analysis at the photon level, we should calculate the cross section for the interesting final state $`2Q\mathrm{\hspace{0.17em}2}\overline{Q}`$, which is by definition an onium-onium scattering. In this case the light quarks in the bottom panel of the Fig. 1 are replaced by a pair of heavy quarks. The production of this particular configuration was for the first time addressed by Szczurek in Ref. , calling the attention that the dipole-dipole approach gives a unique prediction to production of two identical heavy quarks and two identical heavy antiquarks in the final state. The process $`\gamma \gamma 2Q\mathrm{\hspace{0.17em}2}\overline{Q}`$ can only to be produced in next-to-leading order calculation in the collinear approach and/or hadronization process. The cross section for this process is given by,

$`\sigma ^{\mathrm{dd}}(\gamma \gamma 2Q\mathrm{\hspace{0.17em}2}\overline{Q})={\displaystyle |\mathrm{\Psi }_{Q\overline{Q}}(𝒓_1,z_1)|^2|\mathrm{\Psi }_{Q\overline{Q}}(𝒓_2,z_2)|^2\sigma _{\mathrm{dd}}(𝒓_1,𝒓_2,\stackrel{~}{x}_{ab})d^2𝒓_1d^2𝒓_2𝑑z_1𝑑z_2}.`$ (10)
In Fig. 4 are shown the results using the Eq. (10) and considering the saturation model (radius I and II) and the BFKL dipole-dipole cross section for the charm quark case. In this case, $`\overline{r}_0^2=1/m_c^2`$ and we assume $`𝒞=1/9`$ to reproduce the results in the light quark case. Once again, the BFKL result is steeper on energy than the dipole model predictions. Similarly to Ref. , we obtain that at high energies the cross section for $`2c2\overline{c}`$ is about 9 $`\%`$ of that for the single $`c\overline{c}`$ production. It is important to salient that in the calculations of the inclusive $`c\overline{c}`$ production cross section, the contribution of this process should be doubled because each of the quarks can potentially be identified experimentally (See discussion in Ref. ). In next section, we present our results for heavy quark production in ultraperipheral heavy ion collisions, considering both contributions for the inclusive heavy quark production cross section.

## III The heavy quark production in ultraperipheral collisions

Having determined the cross sections for heavy quark production at the photon level, in the following the results for ultraperipheral heavy ions collisions will be presented. We give emphasis on the LHC region, where $`W_{\gamma \gamma }200`$ GeV, and a luminosity $`=4.2\times 10^{26}`$ cm<sup>-2</sup>s<sup>-1</sup> for PbPb collisions ($`\sqrt{s}=5500`$ GeV). Details on the effective two-photon luminosity used here can be found in the previous work in Ref. . As referred above, one has been taken $`m_c=1.3`$ GeV and $`m_b=4.5`$ GeV. Further, we discuss the mass dependence of our results.

Our predictions for the inclusive heavy quark cross section and double $`c\overline{c}`$ pair production at $`\sqrt{s}=5500`$ GeV are presented in Table I. In comparison with the predictions from Ref. , our results for the inclusive cross section are higher by approximately one factor of 1.8. This enhancement is directly associated to the QCD dynamics, which resums higher order diagrams beyond those considered previously in the literature. It is important to salient that this factor is larger than the uncertainty associated with the treatment of the $`\gamma \gamma `$ luminosity, which is estimated to be about 15$`\%`$ . For RHIC energies (not presented in the Table), due to the limited center of mass energy of the photon-photon process, the QCD Pomeron effects are not important, and our results are similar with those presented in Ref. .

In Fig. 5 are presented the results for the energy dependence of the inclusive charm pair production in two-photon ultraperipheral heavy ion collisions. The estimates using the saturation model are shown as the solid and long-dashed lines, corresponding to the choice (I) and (II) for the effective radius, respectively. The dipole BFKL calculation is shown in the dot-dashed curve, where we have considered the Pomeron intercept $`\omega _{IP}=1.47`$. For comparison, the QPM and the single $`c\overline{c}`$ cross sections are also presented. Despite the deviations between the approaches being large at the photon level, the results for the peripheral reaction are quite similar. This feature is expected, since the two-photon luminosity in nuclear collisions strongly suppress high energy contributions in comparison with the behavior present in $`e^+e^{}`$ collisions. Concerning the quark mass dependence of our results, we have checked that if one uses $`m_c=1.2`$ GeV, as done in , the cross section grows approximately 25 $`\%`$. Therefore, for this quark mass choice, we predict an enhancement of a factor two. Here we have used $`m_c=1.3`$ GeV, since with this choice the experimental data are well described in the dipole approach.

In Fig. 6 are shown the results for the energy dependence of the inclusive bottom pair production in two-photon ultraperipheral heavy ion collisions. The solid line represents the estimate obtained with the saturation approach (model I), the dot-dashed line the dipole BFKL result and the dashed line is the QPM prediction. The deviation between the saturation and BFKL approaches reaches up to 10%. Similarly to the charm, the analysis of the bottom quark production in ultraperipheral heavy ion collisions would not allow to discriminate between the distinct kinds of dynamics. However, these effects cannot be disregarded at LHC energies, since our predictions are a factor two larger than those obtained in .

A comment is in order here. In  a detailed analysis of the experimental separation between photoproduction and two-photon interactions was presented. There, the authors have estimated that the two-photon cross sections are at least 1000 smaller than the photoproduction cross section, which makes the experimental separation between the two interactions very hard. Our calculations indicate that the inclusion of the QCD Pomeron effects implies higher cross sections at two-photon level and, consequently, larger cross sections in ultraperipheral collisions. Therefore, the inclusion of these effects implies that the contribution of heavy quark production in two-photon interactions is non-negligible. However, the experimental separation of two-photon process is still remains a challenge.

In Fig. 7 we present separately the results for the double charm pair production in two-photon ultraperipheral heavy ion collisions, obtained from Eq. (10). The solid line represents the estimate obtained with the saturation approach (model I) and the dot-dashed line is the dipole BFKL result. The difference between the predictions of the two approaches is about $`20\%`$ in this particular case. In comparison with inclusive charm pair production, the cross section is two order of magnitude smaller, but still sizeable for an experimental analysis. We salient that the efficiency of flavor tagging is in general very small, which also makes the experimental detection of this process a challenge. However, it important to emphasize that distinctly from the single pair production which is one background of the photoproduction process, the double heavy quark production can be produced in photo-nucleus reactions only if the multiple interactions are considered. Since, in general, these interactions are suppressed by a factor of about 1000 , we have that for this process two photon interactions and photon-nucleus interactions would be of same order.

## IV Summary

We have investigated the heavy quark production in the ultraperipheral heavy ion collisions from two-photon reactions. The total cross sections at the photon level are calculated using the saturation model, describing the dipole-dipole interaction at high energies including in a suitable way the saturation phenomenon. The dipole BFKL approach is also considered as an important test of QCD dynamics and predictions for the photon-photon cross section using this approach are presented. Both calculations are in reasonable agreement with the current experimental measurements. In particular, we have considered a typical final state allowed by the color dipole picture, namely the double heavy quark pair production $`\gamma \gamma 2Q\mathrm{\hspace{0.17em}2}\overline{Q}`$.

The resulting cross sections in ultraperipheral heavy ion collisions are estimated in these distinct dynamical approaches for the QCD Pomeron and presented. Our results demonstrate that the analysis of the heavy quark production in ultraperipheral heavy ion collisions would not allow to discriminate between the distinct dynamics. However, the QCD Pomeron effects cannot be disregarded at LHC energies, since the cross sections are enhanced by a factor two for these effects. Moreover, we have analyzed the double charm pair production and verified that the cross section for this process is sizeable, mainly due to the fact that photoproduction background can only be produced by higher-order multiple interactions. Our results motivate more detailed studies. In particular, the analysis of the mass and transverse momentum distributions could be useful to discriminate between photon-Pomeron and two-photon processes. Work in this direction is in progress.

## ACKNOWLEDGMENTS

The authors thank S. Klein for valuable suggestions and comments. We are particularly grateful to L. Motyka for his careful reading of the manuscript and for his helpful comments. M.V.T.M. thanks the support of the High Energy Physics Phenomenology Group (GFPAE, IF-UFRGS) at the Institute of Physics, Porto Alegre. This work was partially financed by the Brazilian funding agencies CNPq and FAPERGS.