text:
# Onset of Collective Oscillation in Chemical Turbulence under Global Feedback

## I Introduction

Chemical turbulence in oscillatory reaction-diffusion systems can be completely suppressed by means of global feedback ref:battogtokh96 ; ref:battogtokh97 ; ref:kim01 . Theoretically, this fact was found in the complex Ginzburg-Landau equation ref:battogtokh96 ; ref:battogtokh97 , i.e., the normal form of oscillatory reaction-diffusion systems near the supercritical Hopf bifurcation point ref:kuramoto84 . Recent experiments on catalytic CO oxidation on Pt surface demonstrated the same fact, revealing also a variety of wave patterns caused by the effects of global delayed feedback ref:kim01 ; ref:bertram03-1 . A theoretical model for this reaction system reproduced similar behavior ref:kim01 ; ref:bertram03-2 .

In the present paper, we show that yet another transition of universal nature can occur at a certain feedback intensity which is much weaker than the critical intensity associated with the complete suppression of turbulence. The new type of transition is characterized by the emergence of small-amplitude collective oscillation out of the strongly turbulent medium without long-range phase coherence. When the collective oscillation appeared, the system remains strongly turbulent, while the effective damping rate of the uniform mode (i.e., the mean field) shows a change of sign from positive to negative. Thus, the transition is interpreted as a consequence of a complete cancellation of the effective damping of the mean field with the effect of its growth produced by the global feedback.

In Sec. II, we start with the complex Ginzburg-Landau equation (CGL) with global feedback. Then we derive phenomenologically a nonlinear Langevin equation governing the mean field in the form of a noisy Stuart-Landau equation (SL). In order to clarify the nature of the transition of our concern, some numerical results for the one- and two- dimensional CGL will be compared in Secs. III and IV, respectively, with analytical results obtained from the noisy SL. Concluding remarks will be given in the final section.

## II Langevin equation for turbulent CGL as an effective equation

One-dimensional complex Ginzburg-Landau equation with global feedback is given by ref:battogtokh96 ; ref:battogtokh97

$$_tW=W+\left(1+ic_1\right)_x^2W\left(1+ic_2\right)\left|W\right|^2W+\mu \overline{W},$$
(1)
$$\overline{W}\left(t\right)=\frac{1}{L}_0^LW(t,x)𝑑x,$$
(2)
where $`W`$ is a complex field, and $`L`$ is the system size which is supposed to be sufficiently large. The intensity $`\mu `$ of the global feedback is generally a complex number. It is known, however, that a suitable tuning of the delay time in the feedback in the original system can control the phase of this parameter ref:battogtokh96 ; ref:battogtokh97 . For the sake of simplicity, therefore, we shall confine our present analysis to the case of real $`\mu `$, which corresponds to the situation where the delay time in the feedback is fixed at a certain value but the feedback intensity is allowed to vary. A brief comment will be made on the case of complex $`\mu `$ in the final section.

We first consider the system without feedback ($`\mu =0`$), i.e., the usual one-dimensional CGL ref:aranson02 ; ref:shraiman92 . As is well known, uniform oscillations are linearly unstable and turbulence develops when the Benjamin-Feir instability condition

$$1+c_1c_2<0$$
(3)
is satisfied. In what follows, we will fix the parameters $`c_1`$ and $`c_2`$ as $`c_1=2.0`$ and $`c_2=2.0`$ so that the system may stay well within the turbulent regime. We confirmed that under this condition no collective oscillation exists, i.e., $`\overline{W}`$ is randomly fluctuating on a “microscopic” scale around the zero value without perceptible systematic motion. The core of our argument developed below depends little on the choice of parameter values as far as the condition (3) is well satisfied.

It is known that if the turbulence is sufficiently strong, which is actually the case under the above parameter condition, the system exhibits extensive chaos characterized by the property ref:egolf94

$$D_\mathrm{f}L,$$
(4)
where $`D_\mathrm{f}`$ is the Lyapunov dimension of the high-dimensional chaotic attractor describing the turbulence. Extensive chaos implies that the system can be imagined as composed of a large number of cells of equal size such that the fluctuations of some variables associated with the individual cells about their mean value are statistically independent from cell to cell. Thus, the fluctuations of a macro-variable, i.e., a variable given by a simple sum of cell variables over the entire system, are expected to obey the central limit theorem. In particular, in the absence of long-range order, the characteristic amplitude of $`\overline{W}`$ will scale like $`1/\sqrt{L}`$ for large system size. If we represent our continuous oscillatory medium with a long array of $`N`$ oscillators with sufficiently small but fixed separation between neighboring oscillators, which we actually do in numerical simulations to be described below, we expect asymptotically the property

$$\overline{W}=\frac{1}{N}\underset{j=1}{\overset{N}{}}W_j=O\left(\frac{1}{\sqrt{N}}\right),$$
(5)
where $`W_j`$ is the complex amplitude of the $`j`$-th oscillator in the array. The simple nature of the extensive variable $`\overline{W}`$ mentioned above also implies that its time-correlation function defined by

$$C\left(t\right)=\overline{W}^{}(0)\overline{W}(t)\underset{T\mathrm{}}{lim}\frac{1}{T}_0^T\overline{W}^{}(s)\overline{W}(s+t)𝑑s$$
(6)
obeys a simple exponential decay law for large $`t`$, or

$$C(t)\mathrm{exp}\left(\gamma t\right).$$
(7)
The effective damping coefficient $`\gamma `$ is a quantity which could only be determined from the statistical mechanics of turbulent fluctuations which is still far from being established. The above decay law is confirmed from our numerical simulation. Figure 1 shows numerically calculated $`C(t)`$ from which the exponential law with $`\gamma 0.22`$ is confirmed except for some initial transient. In this numerical simulation, and in all numerical simulations which follow, we applied an explicit Euler integration scheme with a constant time step $`\mathrm{\Delta }t0.01`$ and a fixed grid size $`\mathrm{\Delta }x=0.5`$, adopting periodic boundary conditions. The system size $`N`$ used ranges from $`400`$ to $`200000`$.

We now introduce global feedback and study its effects on the dynamics of the mean field. Let the complex amplitude be decomposed into Fourier series as

$$W=\underset{k=\mathrm{}}{\overset{\mathrm{}}{}}\stackrel{~}{W_k}e^{iq_kx},$$
(8)
where $`q_k=2\pi k/L`$ ($`k`$ is an integer), and the Fourier amplitudes are defined by

$$\stackrel{~}{W}_k=\frac{1}{L}_0^LWe^{iq_kx}𝑑x.$$
(9)
The uniform amplitude $`\stackrel{~}{W}_0`$, which is identical with the mean field $`\overline{W}`$ by definition, obeys the equation

$$\dot{\stackrel{~}{W_0}}=\stackrel{~}{W}_0\left(1+ic_2\right)\underset{k_1,k_2}{}\stackrel{~}{W}_{k_1}\stackrel{~}{W}_{k_2}\stackrel{~}{W}_{k_1+k_2}^{}+\mu \stackrel{~}{W}_0.$$
(10)
One may wish to obtain an equation for the mean field in a closed form, which would be a stochastic equation of the nonlinear Langevin type. However, deriving such an equation would be a formidable statistical mechanical problem, so that in the present paper we will content ourselves by simply assuming phenomenologically that the effective exponential decay of the mean field as is seen in Fig. 1 is a result of renormalization of the linear coefficient in Eq. (10) by the nonlinear mode-coupling term. If we can neglect non-Markovian effects, the result of such renormalization will generally take the form of a nonlinear Langevin equation

$$\dot{\overline{W}}=\left[\left(\mu \gamma \right)+i\omega \right]\overline{W}+𝒩\left(\overline{W}\right)+f(t),$$
(11)
where a nonlinear term $`𝒩(\overline{W})`$ with unknown specific form has been included, and $`f(t)`$ represents random force with vanishing mean.

The analysis of Eq. (11) developed in the following section is based on the simplifying assumption that $`f(t)`$ is white Gaussian with the only non-vanishing second moment given by

$$f(t)f^{}(t^{})=4\mathrm{\Gamma }\delta (tt^{}).$$
(12)
The Gaussian nature of the random force seems to hold due to the aforementioned extensive nature of $`\overline{W}`$. The white-noise assumption also seems valid because we are particularly concerned with the situation near the transition point where the characteristic time scale of $`\overline{W}`$ becomes very long. We should also note that the damping coefficient $`\gamma `$ has been assumed to be unchanged when the global feedback is introduced, which is actually the property confirmed by our numerical experiments at least for real $`\mu `$. Our final remark is that in the above argument about Eq. (11) we did not explicitly refer to spatial dimension. Therefore, there seems to be no reason why Eq. (11) should not be applied to systems of two or higher dimensions.

Equation (11) tells that a transition occurs when the global feedback intensity $`\mu `$ becomes equal to the effective damping coefficient $`\gamma `$ of the mean filed. We denote this value of $`\mu `$ as $`\mu _\mathrm{c}`$, or

$$\mu _\mathrm{c}=\gamma .$$
(13)
Figure 2 shows numerically obtained long-time averages of the mean field amplitude, denoted as $`r`$, as a function of the feedback intensity $`\mu `$ over a wide range of $`\mu `$. Although not very clear from these data, there is an indication of transition from vanishing to non-vanishing $`r`$ at some value of $`\mu `$ much smaller than those giving rise to a hysteresis between turbulent and non-turbulent uniform states. In fact, the whole analysis in the rest of this paper is devoted to finding out unambiguous evidence for the existence of a transition through a closer investigation of small-$`\mu `$ region.

In Fig. 3, two space-time plots of the complex amplitude modulus for $`\mu =0.1`$ (weaker feedback) and $`\mu =0.3`$ (stronger feedback) are indistinguishable. In Fig. 4, however, two phase portraits in the complex amplitude plane each obtained for $`\mu =0.1`$ and $`\mu =0.3`$ are contrasted with each other. While the distribution of the representative points of the local oscillators looks almost isotropic for the case of $`\mu =0.1`$, implying the absence of collective oscillations, such symmetry is obviously lost when $`\mu =0.3`$, implying the existence of collective oscillations. Qualitative difference between the two situations is further confirmed from Fig. 5 where a trajectory of the mean field over a long time at each value of $`\mu `$ is displayed. It is clear that when the feedback is weak the mean field is non-oscillatory, simply fluctuating (presumably due to the finite-size effects described above) around the origin, whereas for stronger feedback the same quantity clearly exhibits a closed orbit with some amplitude fluctuation again due to the finite-size effects. Thus, if there is a transition somewhere between these $`\mu `$-values, it is presumably characterized by a noisy Hopf bifurcation.

The nonlinear Langevin equation (11) actually predicts the occurrence of a noisy Hopf bifurcation. From the symmetry of our system and the smallness in amplitude of the collective oscillation near its onset, the nonlinear effects $`𝒩(\overline{W})`$ in Eq. (11) will be dominated by a cubic term. Then the mean field obeys a noisy Stuart-Landau equation

$$\dot{\overline{W}}=\left[\left(\mu \mu _\mathrm{c}\right)+i\omega \right]\overline{W}\left(\eta +i\alpha \right)\left|\overline{W}\right|^2\overline{W}+f(t),$$
(14)
where the coefficients $`\mu _\mathrm{c}`$, $`\omega `$, $`\eta (>0)`$, and $`\alpha `$ are real and depend generally on $`c_1`$ and $`c_2`$. Since the above equation can be handled analytically, it would be interesting to compare some of the results from its analysis with our direct numerical simulations on Eq. (1) and its two-dimensional extension, which are the subjects of Secs. III and IV, respectively.

## III Predicted critical behavior and comparison with numerical results

The nonlinear Langevin equation (14) is equivalent with the Fokker-Planck equation of the following form ref:risken89 ; ref:kampen81 :

$`{\displaystyle \frac{}{t}}P(r,\varphi ,t)`$ $`={\displaystyle \frac{1}{r}}{\displaystyle \frac{}{r}}\left\{\left[\left(\mu \mu _\mathrm{c}\right)r^2+\eta r^4\right]P+\mathrm{\Gamma }r{\displaystyle \frac{}{r}}P\right\}`$
$`+\left[\left(\omega +\alpha r^2\right){\displaystyle \frac{}{\varphi }}P+{\displaystyle \frac{\mathrm{\Gamma }}{r^2}}{\displaystyle \frac{^2}{\varphi ^2}}P\right].`$ (15)
Here $`P(r,\varphi ,t)`$ is the probability density for the amplitude $`r`$ and the phase $`\varphi `$ of the mean field at time $`t`$. The above equation admits a stationary solution independent of $`\varphi `$ given by

$$P_{\mathrm{st}}\left(r\right)=\mathrm{exp}\left[\frac{\eta }{4\mathrm{\Gamma }}\left(r^2\frac{\mu \mu _\mathrm{c}}{\eta }\right)^2\right].$$
(16)
Various moments of $`r`$ defined by

$$r^n=\frac{{\displaystyle _0^{\mathrm{}}}r^nP_{\mathrm{st}}\left(r\right)r𝑑r}{{\displaystyle _0^{\mathrm{}}}P_{\mathrm{st}}\left(r\right)r𝑑r},$$
(17)
and those of the fluctuation $`\delta rrr`$ can be calculated ref:risken65 . In particular, in the limit of weak random force, the mean field amplitude $`r`$ near the transition is found to depend on $`\mu \mu _\mathrm{c}`$ as

$$r=\{\begin{array}{cc}A\left(\mu \mu _\mathrm{c}\right)^{1/2}\hfill & \left(\mu >\mu _\mathrm{c}\right),\hfill \\ 0\hfill & \left(\mu <\mu _\mathrm{c}\right),\hfill \end{array}$$
(18)
where $`A`$ is a constant. Similarly, the mean field fluctuation $`(\delta r)^2`$ is given by

$$(\delta r)^2\left|\mu \mu _\mathrm{c}\right|^1,$$
(19)
which holds for $`\mu \mu _\mathrm{c}`$. It is clear that the critical exponents associated with $`r`$ and $`(\delta r)^2`$ obey the classical law or the mean field theory, reflecting the fact that the transition is caused by the applied mean field and not by the development of local order to a macroscopic scale. Still the transition is in some sense statistical in nature unlike bifurcations in deterministic dynamical systems. This is because the loss of long-range order is solely due to turbulent fluctuations, so that a full theoretical understanding of the transition phenomenon would be impossible without statistical mechanics of chemical turbulence.

Long-time averages of $`r`$ and $`(\delta r)^2`$ were obtained as a function of the feedback intensity from numerical simulation of Eqs. (1) and (2) with $`N=200000`$, and the results are shown in Figs. 6(a) and (b), respectively. By assuming that the long-time averages are identical with ensemble averages, these numerical data were fitted with the theoretical curves given by Eqs. (18) and (19), respectively, where a scale factor is used as the only adjustable parameter. Note that $`\mu _\mathrm{c}`$ is not an adjustable parameter, but is a constant given by Eq. (13).

In obtaining Eq. (18) for the mean field amplitude, we considered the limit of weak random force. Let our argument be generalized to include the dependence of $`r`$ on the noise intensity as well as on $`\mu \mu _\mathrm{c}`$. Because the origin of the random force driving the mean field is the finiteness of the system size, the noise intensity should be inversely proportional to the system size, or

$$\mathrm{\Gamma }=\frac{D}{N},$$
(20)
where $`D`$ is a constant independent of $`N`$. One may thus write the stationary distribution in the form

$$P_{\mathrm{st}}(r)=\mathrm{exp}\left[\frac{N\eta }{4D}\left(r^2\frac{\epsilon }{\eta }\right)^2\right],$$
(21)
where $`\epsilon =\mu \mu _\mathrm{c}`$. Applying the finite-size scaling law developed in Ref. ref:pikovsky99  to the average amplitude of the mean field, we obtain a scaling form

$$r=N^{1/4}F\left(\epsilon N^{1/2}\right),$$
(22)
where $`F`$ is a function (called the scaling function) depending on $`\epsilon `$ and $`N`$ only through $`\epsilon N^{1/2}`$. Equation (22) is a generalization of Eq. (18). Numerically calculated $`r`$ for various $`\epsilon `$ and $`N`$ confirms this scaling law. Figure 7(a) shows the dependence of the long-time average of the mean field amplitude on feedback intensity for some different values of $`N`$. As is seen from Fig. 7(b), all these data come to lie on an identical universal curve after the rescalings of $`r`$ and $`\epsilon `$ by $`N^{1/4}`$ and $`N^{1/2}`$, respectively. In this way, the finite-size scaling law (22) is confirmed, providing unmistakable evidence for a phase transition.

## IV Two-dimensional case

The one-dimensional reaction-diffusion systems which we have numerically studied in Secs. II and III are not very realistic. Especially, surface chemical reactions, such as catalytic CO oxidation on Pt, occur in two dimensions. As mentioned in the foregoing sections, the transition considered here is the mean field type, so that the nature of the transition is expected to be the same as that in one-dimensional systems. We carried out numerical simulations on the two-dimensional complex Ginzburg-Landau equation with global feedback, and obtained a clear indication of transition, although the corresponding value of $`\mu _\mathrm{c}`$ is considerably smaller than that of the one-dimensional case under the same parameter condition. Figure 8 summarizes numerical results for various $`\epsilon `$ and $`N`$, exhibited in a similar manner to Fig. 7(b), i.e., in the form of $`rN^{1/4}`$ vs. $`\epsilon N^{1/2}`$ for different values of $`N`$. These data form almost an identical curve again, which we take as evidence for a transition similar to that in one-dimensional systems.

## V Concluding remarks

Preceding the transition at which the turbulence is completely suppressed and uniform oscillations set in, a different type of transition characterized by the emergence of collective oscillations was shown to exist in the one- and two- dimensional complex Ginzburg-Landau equations with global feedback. The transition is well described phenomenologically with the noisy Stuart-Landau equation governing the mean field. Since the noise there comes from the finite-size effects, the transition becomes infinitely sharp in the limit of infinite system size. The critical exponents of this transition obey the mean field theory because the origin of cooperativity is nothing but the mean field produced by the global feedback. This also gives the reason why spatial dimension one is sufficient for giving rise to the transition. We also confirmed that the two-dimensional complex Ginzburg-Landau equation with global feedback which is more realistic also exhibits a transition of the same type.

Throughout the present paper, our analysis was confined to the case that the feedback intensity $`\mu `$ is real. From our ongoing study, it is being confirmed that a transition of the same nature persists over some range of complex $`\mu `$. The nonlinear Langevin equation (11) also seems to remain valid. Unlike the case of real $`\mu `$, however, the effective damping coefficient $`\gamma `$ appearing in the Langevin equation (11) seems to depend on $`\mu `$, which raises an interesting theoretical problem to be tackled in the future.

Our analysis ref:kawamura06  also suggests that the realistic dynamical model for the catalytic CO oxidation ref:kim01 ; ref:bertram03-2  under experimentally accessible parameter conditions also exhibit a similar transition. We strongly hope for its experimental verification.

###### Acknowledgements.

The authors thank H. Nakao for useful discussions. Numerical computation of this work was carried out with the Computer Facility at Yukawa Institute for Theoretical Physics, Kyoto University.


format: text

text:
# Studies of correlations between 𝐷 and 𝐷̄ mesons in high energy photoproduction

## 1 Introduction

Heavy quark production continues to present itself as a challenge to our understanding of the strong interaction. While Quantum Chromodynamics (QCD) provides a theoretical framework for our understanding and perturbative QCD can be applied to some aspects of heavy quark production, other aspects remain elusive and cannot be described without including a variety of non-perturbative effects. This is especially true for charm production, where perturbative QCD calculations involve large uncertainties and non-perturbative effects play a significant role in modeling physical observables. Until we achieve a fundamental understanding of the strong interaction, accurate models that are able to reproduce properties of the strong interaction—such as heavy quark production—are crucial for our understanding of this fundamental force.

In this paper, we present new results from FOCUS (FNAL-E831) on charm-pair correlations between $`D`$ and $`\overline{D}`$ mesons. Charm-pair correlations have received considerable theoretical attention , and have been studied in both hadroproduction  and photoproduction  experiments. We present our photoproduction results by comparing data distributions to predictions from a recent version of a Monte Carlo based on the Lund Model , which includes non-perturbative effects that have been shown to be important in charm production. We select default settings for charm photoproduction in the Monte Carlo to facilitate comparisons with theoretical predictions and results from other experiments.

## 2 Experimental method

The data for our studies of $`D\overline{D}`$ correlations were recorded by the FOCUS experiment during the 1996–1997 fixed-target run at the Fermi National Accelerator Laboratory. The experiment ran with a photon beam111The photon beam was produced from the bremsstrahlung of secondary electrons and positrons with an endpoint energy of $``$ 300 GeV. The average photon energy for the recorded data was $``$ 180 GeV with a width of $``$ 50 GeV. and a spectrometer that was upgraded from a previous photoproduction experiment, E687 . The FOCUS spectrometer had a target that consisted of four BeO target elements for most of the recorded data222Early in the run a few different targets were used, and less than 5$`\%`$ of the charm-pair data were recorded with Be (instead of BeO) target elements.. A vertex detector, which was located in the target region, had a total of 16 planes of silicon strip detectors. Four of the planes were interleaved with the BeO target elements, and 12 were located downstream of the target. Tracks that were reconstructed in the vertex detector were linked to particle tracks that were found in five multiwire proportional chambers. Particle momenta were determined by measuring the deflection of tracks in two analysis magnets of opposite polarity, and particle identification was accomplished using measurements from three multicell threshold Čerenkov counters, details of which are described elsewhere .

Here we describe, for the first time, the candidate-driven algorithm that was used to collect a large sample of $``$ 7000 pairs of fully reconstructed charmed mesons. The sample consists of pairs of $`D`$ mesons: $`D^+D^{}`$, $`D^+\overline{D}^0`$, $`D^0D^{}`$, and $`D^0\overline{D}^0`$. For this paper, we considered the decay modes $`D^0K^{}\pi ^+`$, $`D^+K^{}\pi ^+\pi ^+`$, $`D^0K^{}\pi ^+\pi ^+\pi ^{}`$, and charged-conjugate modes. The algorithm considered all combinations of two, three, and four charged tracks to find a combination that could be associated with the decay of a single $`D`$ meson, and a second combination of tracks that could be associated with a second $`D`$ decay vertex in the same event. The successful reconstruction of two $`D`$ vertex candidates was followed by the reconstruction of a primary interaction vertex, particle identification cuts, and detachment cuts for the $`D`$ vertices relative to the primary vertex. The goal was to achieve low background levels for each decay mode using a minimum number of cuts.

The first step of the candidate-driven algorithm considers all *pairs* of two-, three-, and four-track combinations in an event. Each combination of tracks represents a possible $`D`$ decay. For each track the algorithm considers all possible combinations of charged *K* or *$`\pi `$* assignments such that the assignments are consistent with the decay of a charged or neutral $`D`$ meson. A particular combination of tracks and the associated particle assignments is referred to as a $`D`$ candidate. The mass of each $`D`$ candidate is calculated using the measured track momenta, and is required to fall within a wide range of 1.6–2.4 $`\text{GeV}/c^2`$. To select events with a $`D`$ and a $`\overline{D}`$, the kaons for the two $`D`$ candidates are required to have opposite charge.

The second step is vertex reconstruction. The goal is to find a pair of $`D`$-decay vertices that can be associated with a primary interaction vertex, and to find all other tracks in the event that can be associated with that primary vertex. This part of the algorithm starts by performing a vertex fit for each $`D`$ candidate. The tracks for each $`D`$ candidate are required to form a vertex with confidence level greater than 1$`\%`$. Pairs of $`D`$ candidates that satisfy the confidence level cut are subjected to two additional vertex cuts. The first cut requires that the momentum vectors of the two candidates intersect with a confidence level greater than 1$`\%`$. The second cut rejects background by rejecting pairs of $`D`$ candidates for which the reconstructed daughter tracks for both $`D`$ candidates form a *single* vertex with confidence level greater than 0.1$`\%`$. This rejects background events in which tracks for both candidates all come from a common vertex. The final phase of the vertex reconstruction treats the two $`D`$ candidates as *seed tracks* to find the primary vertex. Vertex fits are performed by including the two seed tracks as well as combinations of all other tracks in the event. As many tracks as possible are added to the primary vertex as long as the confidence level is greater than 1$`\%`$.

Pairs of $`D`$ candidates that survive the vertex reconstruction are subjected to particle-identification cuts, which are based on measurements from three multicell threshold Čerenkov counters. The Čerenkov algorithm  calculates four likelihoods that correspond to the four hypotheses (electron, pion, kaon, proton) that are considered for each charged track. The algorithm produces a $`\chi ^2`$-like variable $`W_i=2\mathrm{ln}(\mathrm{likelihood})`$, where $`i`$ is the index used to represent each hypothesis. For the kaon in each $`D`$ candidate, we require that the kaon hypothesis is favored over the pion hypothesis by more than a factor of exp(0.5) by requiring $`W_\pi W_K>1.0`$. For the pions in each $`D`$ candidate we apply a pion consistency cut, which requires that no particle hypothesis is favored over the pion hypothesis with a $`\mathrm{\Delta }W=W_\pi W_{\text{min}}`$ greater than $`5`$, where $`W_{\text{min}}`$ is the $`W_i`$ with the smallest value.

After applying particle-identification cuts, we impose cuts based on the significance of detachment ($`\mathrm{}/\sigma _{\mathrm{}}`$) between each $`D`$ candidate and the primary vertex. We calculate $`\mathrm{}/\sigma _{\mathrm{}}`$ by using the measured value of $`\mathrm{}`$, the distance between the $`D`$ decay vertex and the primary vertex, and dividing by the associated error $`\sigma _{\mathrm{}}`$. The cuts for $`\mathrm{}/\sigma _{\mathrm{}}`$ range from $`\mathrm{}/\sigma _{\mathrm{}}>1`$ to $`\mathrm{}/\sigma _{\mathrm{}}>4`$ depending on the decay mode, whether the $`D`$-decay vertex is located between target elements (for which background levels are low) or in target material, and whether a $`D`$ candidate can be associated with a $`D^{}`$ decay.

Figure 1a shows the $`D\overline{D}`$ signal that we obtain after all of the aforementioned cuts have been applied to the FOCUS data. Figure 1a shows the normalized $`D`$ invariant mass333The normalized mass, $`M_n(D)=\mathrm{\Delta }M/\sigma _M`$, is defined as the difference between the reconstructed mass and the central value of the $`D^+`$ or $`D^0`$ mass distribution divided by the reconstructed-mass error $`\sigma _M`$, which is calculated for each $`D`$ candidate. $`M_n(D)`$ opposite the normalized $`\overline{D}`$ invariant mass $`M_n(\overline{D})`$. Figure 1b shows a Gaussian fit to $`M_n(\overline{D})`$ over a linear background after applying a background subtraction procedure that is used to determine the number of pairs of charmed $`D`$ mesons in the FOCUS data. The procedure consists of performing a sideband subtraction and fit for one normalized mass distribution by selecting entries in the signal and sideband regions of the other normalized mass distribution. In Fig. 1b we plot $`M_n(\overline{D})`$ by assigning unit weight to $`D`$ candidates with a reconstructed mass in the signal region ($`\pm 2\sigma `$ about the central value of the $`D^+`$ or $`D^0`$ mass of the candidate), and a weight of $`^1/_2`$ to candidates with mass in the two 4–8$`\sigma `$ sideband regions444An equivalent approach to determine the $`D\overline{D}`$ yield is a fit to $`M_n(D)`$ after selecting signal and sideband regions for $`M_n(\overline{D})`$. Using this approach we obtain a $`D\overline{D}`$ yield of $`7126\pm 120`$ , which is consistent with the yield mentioned in the text.. The $`D\overline{D}`$ yield that we obtain from our fit is $`7064\pm 119`$ (statistical error).

In addition to our study of correlations between pairs of fully reconstructed $`D`$ mesons, we study correlations between two $`D`$ mesons where one $`D`$ is fully reconstructed and the other is kinematically tagged by a slow pion coming from the decay $`D^+\pi ^+D^0`$. In these decays, the $`D^0`$ need not be reconstructed, and therefore we refer to this sample of charmed $`D`$ mesons as partially reconstructed charm pairs555The partially reconstructed sample consists of $`D^+D^{}`$, $`D^+\overline{D}^0`$, $`D^0D^{}`$, and $`D^+D^{}`$ pairs. . The reason for including this sample in our studies of correlations is that charm-pair correlations can be studied over a larger kinematic range compared to the fully reconstructed sample.

For partially reconstructed charm pair events  we begin by considering all two-, three-, and four-track combinations for the fully reconstructed $`D`$ (recoil $`D`$) in an event. We consider the decay modes $`D^0K^{}\pi ^+`$, $`D^+K^{}\pi ^+\pi ^+`$, $`D^0K^{}\pi ^+\pi ^+\pi ^{}`$, and charged-conjugate modes. A candidate-driven algorithm uses the recoil $`D`$ candidates to find the primary vertex, requiring the vertex confidence level to be greater than 1$`\%`$. The same Čerenkov particle identification criteria used for the fully reconstructed charm-pair sample (see above) are applied to the recoil $`D`$ candidates. However, a more restrictive detachment cut of $`\mathrm{}/\sigma _{\mathrm{}}>5`$ is applied to all three decay modes in the partially reconstructed charm-pair sample. Figure 1c shows the invariant mass distribution, which includes all three decay modes, with a total of $`\mathrm{782\hspace{0.17em}630}\pm 1600`$ candidates satisfying the selection criteria.

The next step in the analysis treats each track that is assigned to the primary vertex (excluding the recoil $`D`$) as a slow-pion candidate from the decay $`D^+\pi ^+D^0`$. The momentum of the track is multiplied by 13.8 to approximate the momentum of the $`D^+`$666Due to the low $`Q`$ value of the $`D^{}`$ decay, the momentum of the soft pion approximates the momentum of the $`D^{}`$ when multiplied by the inverse of its energy fraction, which is $``$ 13.8. . If the charge of the slow pion is the same as the charge of the kaon from the recoil $`D`$, then the combination of the slow pion and recoil $`D`$ is designated as a *right-sign* combination. Otherwise, it is a *wrong-sign* combination. This assignment of right- and wrong-sign combinations is used for background subtraction.

A double subtraction method is used to reduce backgrounds. First, to handle non-charm background, a sideband subtraction is applied to recoil $`D`$ candidates. A Gaussian fit is applied to the invariant mass distribution for each of the three decay channels. Entries in the 4–8$`\sigma `$ sideband regions are subtracted from those in the $`\pm 2\sigma `$ peak region by using a weight factor of $`^1/_2`$. Second, the assignment of right- and wrong-sign combinations is used to subtract wrong-sign background from right-sign combinations. To avoid distortion of the wrong-sign background we *exclude* all slow-pion candidates that can be associated with a $`D^{}`$ decay involving the recoil $`D`$. This *anti-$`D^{}`$* cut is imposed by excluding combinations of slow pions and recoil $`D`$-mesons that have a mass difference, $`m(D^{})m(D)`$, in the range 0.142–0.149 $`\text{GeV}/c^2`$. To further enhance the selection procedure, a maximum cut of 4 (GeV/$`c`$)2 is applied to $`\mathrm{\Delta }_t^2=(p_x^{(r)}+13.8p_x^{(\pi )})^2+(p_y^{(r)}+13.8p_y^{(\pi )})^2`$, where $`p_x^{(r)},p_y^{(r)}`$ and $`p_x^{(\pi )},p_y^{(\pi )}`$ are transverse momentum components of the recoil $`D`$ and slow pion, respectively. This cut enhances the selection of signal since genuine events balance $`\mathrm{\Delta }_t^2`$ (see reference  for more details). This is shown in Fig. 1d, which shows a prominent excess of right-sign combinations close to $`\mathrm{\Delta }_t^2=0`$ compared to the wrong-sign background. After applying the double subtraction and the $`\mathrm{\Delta }_t^2`$ cut, we obtain a sample of $`\mathrm{75\hspace{0.17em}160}\pm 1040`$ partially-reconstructed charm pairs.

## 3 $`D\overline{D}`$ correlations

For our study of correlations between pairs of fully reconstructed $`D`$ mesons, we compare FOCUS data to predictions from a Monte Carlo based on the Lund Model. The Monte Carlo consists of a Pythia 6.203  generator with default settings, and detector simulation algorithms for the FOCUS apparatus. The Monte Carlo generator produces charm events using a tree-level photon-gluon fusion process applied to beam photons and target nucleons. We use default options for charm photoproduction in the generator (instead of using a Monte Carlo tuned to match our data) to facilitate comparisons with theoretical predictions and results from other experiments. In this paper, we also compare our results to previously published charm photoproduction results from experiment E687 .

To improve comparisons between data and model predictions based on photon-gluon fusion, we eliminate our lowest multiplicity events by requiring a minimum number of particles assigned to the primary interaction vertex. We define $`N_{\text{primary}}`$ as the number of particles assigned to the primary vertex. With this definition, $`N_{\text{primary}}`$ has a minimum value of two since it includes the $`D`$ and $`\overline{D}`$ mesons (each charm meson counts as a single particle) in addition to charged tracks assigned to the primary vertex. To eliminate our lowest multiplicity events we require an $`N_{\text{primary}}>2`$ cut. The cut eliminates features observed in data that are not present in Pythia 6.203. This is illustrated in Fig. 2a, which shows the background-subtracted777The background subtraction procedure assigns unit weight to $`D\overline{D}`$ candidates in the signal region in Fig. 1a ($`\pm 2\sigma `$ about the center of the distribution), a weight of $`^1/_2`$ to candidates in the single $`D`$ and single $`\overline{D}`$ sidebands (four regions defined as $`\pm 2\sigma `$ about the $`D`$ axis and $`\pm `$ 4–8$`\sigma `$ about the $`\overline{D}`$ axis, and $`\pm 2\sigma `$ about the $`\overline{D}`$ axis and $`\pm `$ 4–8$`\sigma `$ about the $`D`$ axis), and a weight of $`+^1/_4`$ to candidates in the four regions where both the $`D`$ and $`\overline{D}`$ candidates are 4–8$`\sigma `$ away from the center of the distribution. The weight factor of $`+^1/_4`$ accounts for the over-subtraction of the single-$`D`$ and single-$`\overline{D}`$ backgrounds and the subtraction of random combinatoric background. invariant $`D\overline{D}`$ mass for mass combinations with a net charge of zero ($`D^+D^{}`$ and $`D^0\overline{D}^0`$) for FOCUS data, and for Pythia events that have passed through a software simulation of the FOCUS detector and have survived the event selection procedure described earlier in this paper. The mass distribution has an enhancement near threshold that is not present in Pythia. This enhancement is evident for events with $`N_{\text{primary}}=2`$, especially when we apply additional cuts that remove events with energy deposited in electromagnetic calorimeters (see inset in Fig. 2a). The enhancement seems to arise from the diffractive production of $`\psi (3770)`$ decaying to $`D\overline{D}`$, and will be the subject of a future paper (additional information can be found in conference proceedings ). Another significant difference between data and Pythia is shown in Figs. 2a and 2b. This is the excess of $`N_{\text{primary}}=2`$ events in data compared to Pythia, some of which can be attributed to the production of $`\psi (3770)`$. By eliminating $`N_{\text{primary}}=2`$ events, we get fairly good agreement for the $`N_{\text{primary}}`$ distribution in Fig. 2b, which shows the histogram for Pythia (solid line) normalized to the number of $`D\overline{D}`$ pairs in the data (data points with error bars) with $`N_{\text{primary}}>2`$. By eliminating the $`N_{\text{primary}}=2`$ bin the agreement between data and Pythia is significantly improved (a slight excess of events with $`N_{\text{primary}}=3`$ persists in the data).

Previous studies  of charm-pair correlations have presented distributions for $`p_t^2(D\overline{D})`$, the transverse momentum squared of the $`D\overline{D}`$ pair, and $`\mathrm{\Delta }\varphi `$, the azimuthal angle between the $`D`$ and $`\overline{D}`$ momentum vectors in the plane transverse to the beam direction. These distributions are significant, since $`p_t^2(D\overline{D})=0`$ and $`\mathrm{\Delta }\varphi =\pi `$ radians in leading-order QCD, where the charm-quark pair is produced back-to-back. In QCD these distributions are broadened by NLO corrections and non-perturbative effects, as illustrated in references  and . Photoproduction results from E687  have been compared to results from NLO calculations  and Pythia version 5.6 . The E687 comparisons between data and Pythia 5.6 are reproduced in Figs. 2c and 2d, but with a different normalization to match FOCUS data (shown as open circles with error bars). The figures show good agreement between FOCUS and E687 data, and a significant discrepancy between data and Pythia 5.6.

Agreement between FOCUS data and the more recent Pythia 6.203 is significantly better, but minor discrepancies persist. Figure 3 shows comparisons for $`\mathrm{\Delta }\varphi `$, $`p_t^2(D\overline{D})`$, rapidity difference defined as $`\mathrm{\Delta }y`$ = $`y_Dy_{\overline{D}}`$, and invariant $`D\overline{D}`$ mass, $`M(D\overline{D})`$. FOCUS data are plotted as data points with error bars. Pythia parent distributions (dashed lines) are shown without acceptance or resolution effects, so that parent distributions can be compared to the distributions that are obtained for Monte Carlo events that have survived detector simulation, event selection and analysis cuts (solid histograms).

Figure 3a shows good agreement for $`\mathrm{\Delta }\varphi `$. There is an enhancement in the first $`\mathrm{\Delta }\varphi `$ bin, which is not present in Pythia and may suggest the presence of an additional production mechanism. There is good agreement for $`p_t^2(D\overline{D})`$ in Fig. 3b, except that the data tend to have slightly larger values of $`p_t^2(D\overline{D})`$. Compared to Pythia 5.6, the agreement between data and Pythia 6.203 for $`\mathrm{\Delta }\varphi `$ and $`p_t^2(D\overline{D})`$ is significantly better. Some of the improvement can be attributed to a larger value for the intrinsic transverse momentum of the incoming partons, referred to as the $`k_T`$ kick888A value of $`<k_T^2>`$ = $`(1\text{GeV}/c)^2`$ was introduced with Pythia version 6.135, while previous versions had a value of $`<k_T^2>`$ = $`(0.44\text{GeV}/c)^2`$, but a number of other Pythia modifications that affect these distributions have also occurred over time. Figure 3c shows fairly good agreement for $`\mathrm{\Delta }y`$999The agreement between data and Pythia improves slightly for $`D`$ mesons with larger values of $`\mathrm{}/\sigma _{\mathrm{}}`$, however a more restrictive $`\mathrm{}/\sigma _{\mathrm{}}`$ cut also reduces the number of charm-pair events that are available for correlation studies., but also shows significant acceptance losses for $`|\mathrm{\Delta }y|>1`$ (acceptance losses are less severe in the partially reconstructed charm-pair sample). Acceptance losses are also significant for large values of $`M(D\overline{D})`$ in Fig. 3d, but here there is a discrepancy between data and Pythia for smaller values of $`M(D\overline{D})`$ where the acceptance is good.

Figure 4 shows results for the partially reconstructed charm-pair sample, comparing data (asterisks with error bars) to Pythia 6.203 (solid lines). For $`\mathrm{\Delta }\varphi `$ (see Fig. 4a) we also include a comparison to the distribution that we obtain for fully reconstructed charm pairs after accounting for resolution broadening effects101010The $`\mathrm{\Delta }\varphi `$ distribution for the fully-reconstructed sample is obtained by taking the momentum vector of the $`D`$ or $`\overline{D}`$ in an event and treating it as the momentum of a $`D^{}`$ that decays isotropically to a $`D^0`$ and a pion. The pion momentum vector is then used to determine $`\mathrm{\Delta }\varphi `$ as is done in the analysis of partially reconstructed charm pair events.. This shows that the two samples are in agreement, and that the enhancement that we observe in the first $`\mathrm{\Delta }\varphi `$ bin for fully reconstructed charm pairs (see Fig. 3a) disappears due to resolution broadening and selection cuts applied to the partially reconstructed charm-pair sample. Figures 4a and 4b are both affected by resolution broadening (the effects are reproduced by our Monte Carlo), and the agreement between data and Pythia 6.203 is good. As before, the data tend to have slightly larger values of $`p_t^2(D\overline{D})`$.

In Figs. 4c and 4d we show results for $`\mathrm{\Delta }y`$ and $`M(D\overline{D})`$, and include Pythia parent distributions (dashed lines) to show how acceptance losses in this sample compare to acceptance losses in the fully-reconstructed sample (see Fig. 3). The partially reconstructed charm pairs are less affected by acceptance losses, and thus extend the kinematic range of our correlation studies. The $`\mathrm{\Delta }y`$ distributions in Fig. 4c show good agreement, while the $`M(D\overline{D})`$ distributions in Fig. 4d exhibit a mismatch between data and Pythia that is similar to the mismatch that is observed in Fig. 3d.

## 4 Conclusions

We have extracted two large samples of photoproduced charm-pair events for studies of correlations between $`D`$ and $`\overline{D}`$ mesons. The first sample consists of more than 7000 fully reconstructed $`D\overline{D}`$ pairs. The second sample consists of over $`\mathrm{75\hspace{0.17em}000}`$ partially reconstructed charm pairs, where one $`D`$ meson is fully reconstructed and the other is tagged by a slow pion coming from a $`D^{}`$ decay. For the fully reconstructed sample we impose an $`N_{\text{primary}}>2`$ cut to eliminate our lowest multiplicity events, while the partially reconstructed sample has an implicit cut of $`N_{\text{primary}}>2`$ due to the presence of the slow pion. The significance of the $`N_{\text{primary}}`$ cut is that it improves our comparisons to model predictions based on photon-gluon fusion by eliminating low multiplicity events in which we observe the production of $`\psi (3770)`$ decaying to $`D\overline{D}`$ pairs. The $`\psi (3770)`$ events, which are not included in Pythia, appear to be produced diffractively, and will be the subject of a future paper.

The FOCUS results on charm-pair correlations presented in this paper are in good agreement with previous measurements from experiment E687, which displayed significant discrepancies compared to an older version of Pythia (version 5.6). Comparisons of FOCUS data to a more recent version of Pythia (version 6.203) are significantly better, due to changes in parameters that affect the modeling of photon-gluon fusion. One notable change that improves the agreement with data is that the intrinsic transverse momentum ($`k_T`$) of incoming partons was increased from $`<k_T^2>`$ = $`(0.44\text{GeV}/c)^2`$ to $`<k_T^2>`$ = $`(1\text{GeV}/c)^2`$. Although minor discrepancies persist when FOCUS data are compared to Pythia, the modeling of heavy quark photoproduction is fairly good for correlations between $`D`$ and $`\overline{D}`$ mesons.

## 5 Acknowledgments

We wish to acknowledge the assistance of the staffs of Fermi National Accelerator Laboratory, the INFN of Italy, and the physics departments of the collaborating institutions. This research was supported in part by the U. S. National Science Foundation, the U. S. Department of Energy, the Italian Istituto Nazionale di Fisica Nucleare and Ministero dell’Istruzione dell’Università e della Ricerca, the Brazilian Conselho Nacional de Desenvolvimento Científico e Tecnológico, CONACyT-México, and the Korea Research Foundation of the Korean Ministry of Education.


text:
# 1. Introduction

## 1. Introduction

Nanographite systems, where graphene sheets of the orders of the nanometer size are stacked, show novel magnetic properties, such as, spin-glass like behaviors , and the change of ESR line widths while gas adsorptions . It has been found  that magnetic moments decrease with the decrease of the interlayer distance while water molecules are attached physically. Such the phenomena indicate the tunable magnetism in nanometer size systems. Recently, unique magnetisms in carbon-based materials have been observed , and thus theoretical investigations are urged in order to resolve the mechanisms.

This paper reviews recent theoretical works on electronic properties in nanographite materials, and is organized as follows. (1) Antiferromagnetism in the stacked nanographite is investigated with using the Hubbard-type model taking into account of the hopping interactions of $`\pi `$-electrons and the strong onsite repulsions at carbon atoms . We point out the importance of the open shell electronic structure as an origin of the decreasing magnetic moment with the decrease of the inter-layer distance . (2) Possible charge-separated states in nanographite ribbons are discussed in terms of the extended Hubbard model with nearest-neighbor interactions . Such the charge-polarized state could be observed, when a static electric field is applied in the direction parallel with the graphene plane for example. (3) Superperiodic patterns with a long distance in a nanographene sheet observed by STM are discussed in terms of the interference of electronic wave functions. The period and the amplitude of the oscillations decrease spatially in one direction. We explain the superperiodic patterns with a static linear potential in the $`𝒌𝒑`$ model .

## 2. Tuning magnetism in stacked nanographite

First, we report magnetic properties for the A-B stacked hexagonal nanographite shown in Fig. 1 (a) . The first and second layers are displayed by the thick and thin lines, respectively. In each layer, the nearest neighbor hopping $`t`$ is considered. Each layer has closed shell electron systems when the layers do not interact mutually, because the number of electrons is even and equal to the number of sites. The interlayer hopping $`t_1`$ is assigned at the sites with closed circles. The model is solved with the unrestricted Hartree-Fock approximation, and antiferromagnetic solutions are obtained. Figure 1 (b) shows the absolute value of the total magnetic moment per layer as functions of $`t_1`$ and $`U`$. As increasing $`U`$, the magnitude of the magnetization increases. The magnetic moment is zero at the smaller $`t_1`$ region for $`U=1.9t`$ (open squares), $`2.0t`$ (closed circles), and $`2.1t`$ (open circles). The magnetic moment is zero only at $`t_1=0`$ for $`U=2.2t`$ (closed triangles) and $`2.3t`$ (open triangles). We can understand the parabolic curves as a change due to the Heisenberg coupling proportional to $`t_1^2/U`$.

Next, we consider systems with open shell electronic structures when a nanographene layer is isolated . One case is the effects of additional charges coming from functional side groups with introducing a site potential $`E_s`$ at edge sites. We take $`E_s=2t`$, and one additional electron per layer is taken account. Figure 2 displays the absolute values of total magnetic moment per layer. In Fig. 2 (a), the site potentials locate at the site A in the first layer \[Fig. 1 (a)\], and at its symmetrically equivalent site in the second layer. In Fig. 2 (b), the site potential exists at the site B. The total magnetization is a decreasing function of $`t_1`$ in both figures. The decrease is faster in Fig. 2 (b) than in Fig. 2 (a). The site B is neighboring to the site with the interaction $`t_1`$, and thus the localized character of the magnetic moment can be affected easily in this case. The decease of observed magnetization by the magnitude $`3040`$% with the water molecule attachment  may correspond to the case of Fig. 2 (b). The other origin, i.e., the geometrical effects can give rise to the decrease of the magnetic moment with increasing interlayer hopping interactions, too. Therefore, we conclude that the open shell nature of electronic states in a single graphene layer plays an important role in controlling the magnetism in nanographite materials.

## 3. Charge- and spin-separated states in nanographite ribbons

In section 2, we have discussed the magnetic properties with using the onsite interactions $`U`$ only. In this section, we will report a possible charge-separated state in nanographite ribbons with zigzag edges by introducing the nearest-neighbor Coulomb repulsion term $`V`$ as well as $`U`$ . The model is solved with a mean field method with the finite space geometry shown in Fig. 3. Charge- and spin-polarized solutions are obtained depending on the interaction strengths. For example, Fig. 3 (a) shows the charge density distribution of the ferroelectric charge separated (CS) state with $`U=0.3t`$ and $`V=0.4t`$. This state appears and becomes stable when the effects of $`V`$ overcom those of $`U`$. On the other hand, Fig. 3 (b) shows the spin density profile of the localized magnetic (LM) state for $`U=t`$ and $`V=0`$. This state is stable when $`V`$ is weak. The CS (LM) state has no spin (charge) density at every site. In the CS state, the upper (lower) zigzag edge is charged positively (negatively). This distribution pattern is quite similar to that in the LM state. Such the static charge polarization (magnetization) is explained by the interplay of the Fermi instability of the flat bands due to $`V`$ ($`U`$) and the localized edge states. The signs of the charge (spin) densities at neighboring sites are opposite, reflecting the bipartite nature.

Figure 4 shows the phase diagram in the parameters $`U`$ and $`V`$, representing the stability between the CS and LM states. We use the geometry with $`4\times 40`$ carbon atoms. Above the phase boundary, the CS state has the lower energy, and the LM state becomes stable below the boundary. The phase boundary indicates the first order phase transition. At $`U=0`$, the phase boundary rises up with the infinite gradient. This is a signature of the localized edge states, differently from the graphite sheet. The inset shows the phase diagram in reduced scales. In the strong correlation limit, the phase boundary approaches to an analytical phase boundary line $`V=[N/(3N1)]U`$ between the CDW and SDW states, which is obtained by equalizing their Coulomb energies. Here, $`N`$ is the number of zigzag lines of the ribbon. The CDW and SDW solutions are crossovered from those of CS and LM for the weak interactions. Such the charge-separated states could be observed when static electric field is applied in the nanographene plane for example, and the presence of the state will give effects on dielectric properties.

## 4. Electronic wave interference patterns

In Fig. 5, an STM image of the graphene sheet with a necktie shape is shown . The observation has been done with the following condition: bias voltage $`V=200`$ mV and current $`I=0.7`$ nA. The distance between the graphene necktie and the substrate is over 0.8 nm, suggesting that it consists of a stacking of two graphene layers, which interact weekly with the HOPG substrate. Interestingly, the period and the amplitude of the oscillations decrease from the top to the bottom along the graphene necktie. The oscillation period is one order of magnitude larger than that of the Moiré pattern due to stacking, which has been reported elsewhere, and therefore this possibility may be weak. We can assume effects of long-distance periodic-structural deformations in the graphene surface or interference effects of electronic wave functions. We have also observed that the oscillation period becomes longer by placing a nanographene flake on the graphene necktie. The oscillations period seems to be double in the upper region of the necktie after addition of one flake. The oscillation below the flake seems to be only slightly modified by the flake. Such effect on the oscillations cannot be explained by some structural modulations. Therefore, the oscillation patterns could be the effect of interference of the electronic wave functions in the graphene surface.

In order to analyze the interference patterns, we give comparison with the calculation of the model for the graphene plane. Here, we use the continuum $`𝒌𝒑`$ model . The electron density is calculated with including a static potential which has a functional form of the linear decrease from top to the bottom along the surface of the graphene necktie. The peak positions of the electron density in the long direction of the graphene necktie of Fig. 5 are plotted in Fig. 6, and comparison with the theoretical results is given. The decrease of the oscillation period fairly agrees with the experiments. However, the slight decrease of the experimental corrugation cannot be reproduced by the theoretical result because of neglecting an effect of a tip-apex shape of STM on the observed corrugatrion amplitude. The fitting gives the parameter of the potential gradient $`6.49\times 10^3`$ eV/nm. The total potential variation over the distance 200 nm becomes 1.3 eV. Such magnitude of the potential change would survive under thermal lattice fluctuations and can really exist in experiments. The present result by no means implies that the wave functions observed with superperiodic amplitudes are of the electrons which have energy levels of the graphene plane.

## 5. Summary

First, antiferromagnetism in the stacked nanographite has been investigated with the Hubbard-type model. The A-B stacking is favorable for the hexagonal nanographite with zigzag edges, in order that magnetism appears. We have also found that the open shell electronic structures can be origins of the decreasing magnetic moment with adsorption of molecules.

Next, possible charge-separated states have been considered using the extended Hubbard model with nearest-neighbor interactions. The charge-polarized state could appear, when a static electric field is applied in the graphene plane for example.

Finally, we have characterized theoretically the superperiodic patterns in a nanographene sheet observed by STM. We have adopted the $`𝒌𝒑`$ model for the description of the electronic structures of the graphite. The calculated electron density has the property that the oscillation period decreases spatially while the amplitude remains constant. The magnitude of the static potential seems reasonable. It turned out that the long distance oscillations come from electrons with the band structures of the two dimensional graphene sheet.

Figure Captions

Fig. 1. (a) A-B stacked hexagonal nanographite with zigzag edges. (b) The absolute magnitude of the total magnetic moment per layer as a function of $`t_1`$. The onsite interaction is varied within $`1.8t`$ (closed squares) $`U2.3t`$ (open triangles). The interval of $`U`$ between the series of the plots is $`\mathrm{\Delta }U=0.1t`$.

Fig. 2. The absolute magnitude of the total magnetic moment per layer as a function of $`t_1`$ for the system with a site potential $`E_s=2t`$, (a) at the site A and (b) at the site B. The site positions are displayed in Fig. 1 (a). In (a), the onsite interaction is varied within $`0.6t`$ (closed squares) $`U1.8t`$ (closed triangles) with the interval $`\mathrm{\Delta }U=0.3t`$. In (b), it is varied within $`1.0t`$ (closed squares) $`U2.0t`$ (closed triangles) with the interval $`\mathrm{\Delta }U=0.25t`$.

Fig. 3. (a) Charge density distribution of the charge-separated (CS) state, and (b) the spin density distribution of the localized magnetic (LM) state on a zigzag ribbon with $`4\times 20`$ carbon atoms, where $``$ means positive charge (spin) densities, and $``$ indicates negative charge (spin) densities. The radius of each circle denotes the magnitude of the density.

Fig. 4. The phase diagram in the parameter space of $`U`$ and $`V`$ for a zigzag ribbon with $`4\times 40`$ carbon atoms. The solid curve is the boundary between the CS and LM states, interpolating the numerical data ($``$). The inset shows the phase diagram in reduced scales. The dashed line denotes the phase boundary in the strong correlation limit.

Fig. 5. STM image of the superperiodic pattern observed on a necktie shaped graphene plate on HOPG substrate.

Fig. 6. Comparison for the electron wave patterns by STM and the $`𝒌𝒑`$ model. Experimental peak positions along the perpendicular direction of Fig. 5 are plotted by diamonds. The results of the fitting by the long distance envelope functional form derived from the $`𝒌𝒑`$ model are shown by squares.


text:
# Probing human response times

## 1 Introduction

There can be little doubt that the World-Wide-Web (WWW) and Internet e-mail provides two of the most efficient methods for retrieving and distributing information. As such, it carries an enormous potential with respect to sale and marketing of all kinds of products. A rather troublesome feature from a research perspective of many “old-age” communication channels, such as newspapers, radio and TV, interpersonal contacts and so forth, is that the diffusion of information through these channels is in general very difficult and/or time consuming to probe. This is not so with the WWW and e-mail. The fact that these communication channels are computer-based and access in principle unrestricted provides a rather unique opportunity to study in real time how fast individuals respond to a new piece of information. It is generally believed that the rate of diffusion of information in a given population depends on a number of socio-economical factors such as education, cultural status, exposure to mass media and interpersonal channels etc. . The evidence presented here propose otherwise.

Along a complementary line of research the WWW has provide a similarly unique opportunity to study fast evolving (sociological or computer) networks . Most studies of the WWW have until now focused on the more easily accessed statistical properties such as the connectivity of the WWW (the distributions of outgoing and incoming links etc.) and a consensus that the WWW is scale-free has been established. As the WWW or any other sociological network is constantly evolving any such investigation can only deliver a “snap-shot” picture of the network. The assumption is then that the dynamical features of the network is reasonably stable over the time-scales considered and that the results are representative over time. However, in order to correctly describe networks between interacting humans one must estimate characteristics specific to the nature of the “nodes”, e.g., psychological traits of humans, as one may otherwise be mislead when generalizing from a “snap-shot” analysis. In fact, little is known about the internaut population’s response to some external or internal event. With respect to commercial exploitation of the WWW it is obviously the response of the internaut population to some new piece of information which is of prime interest and not the topology or connectivity of the underlying network itself.

The purpose of the present paper is to provide for a characterization of the psychological/sociological traits of interacting humans in absence of formal deadlines. Using empirical data of response time distributions of internaut individuals/populations when exposed to a new piece of information, a surprisingly robust “law” for human response times will be put forward.

## 2 Response to an Internet Interview: Experiment I

The Nasdaq crash culminating on Friday the 14th of April 2000 caught many people with surprise and shook the stock market quite forcefully. As always, many different reasons for the crash were given ranging from the anti-trust case against MicroSoft to an “irrational exuberance” of the participants on the stock market. The author also had bid for the cause , which was made public on Monday the 17th of April 2000 on the Los Alamos preprint server. As a result, a forty minute interview with the author called “The World (Not) According to GARCH” was published on Friday the 26th of May 2000 on a “radio website” (www.wallstreetuncut.com) In addition, the URL to the author’s papers was announced making it clear that work on stock market crashes in general and the recent Nasdaq crash in particular could be found using the posted URL. As the main subject of the interview was a somewhat technical discussion about whether stock market crashes could be predicted, it seems reasonable to assume that the population in a socio-economic context was relatively homogeneous. The experiment proceeded as follows. The number of downloads of papers from the authors homepage as a function of time (days) from the appearance of the interview was recorded. In figure 1, the cumulative number of downloads $`N\left(t\right)`$ as a function of time $`t`$ after the appearance of the interview is shown on a semi-logarithmic scale. The data is over one and half decade surprisingly well-captured by the relation

$$N\left(t\right)a\mathrm{ln}\left(t/b\right)+ctn\left(t\right)=\frac{dN\left(t\right)}{dt}\frac{1}{t}+c.$$
(1)
Here $`k`$ represents a constant background rate. After approximately $`60`$ days, the data breaks away from the fitted line. The reason is the appearance of the author’s URL on the Social Science Research Network server (www.ssrn.com) causing a second advertisement. Hence, the experiment became influenced by an additional distribution channel and was consequently halted after 69 days.

## 3 Response times in Email Interchange: Experiment II

In a recent preprint , the temporal dynamics of an e-mail network has been investigated by J.P. Eckmann, E. Moses and D. Sergi. Specifically, the time period between an e-mail message and its reply was recorded. The data set contains 3188 users interchanging 309129 messages and was obtained from the log-files from one of the main mail servers of an university. In one sense, this experiment is the “cleanest” as only a single distribution channel can exist. However, as each email message contains different information it is also the “dirtiest” with respect to the information distributed. Furthermore, it seems reasonable to assume that the population in a socio-economic context is quite heterogeneous opposed to that of the first experiment. In figure 2, we see that the cumulative distribution of time periods between an e-mail messages and its reply $`N(t)`$ can be modeled in terms of a response time $`t`$ of the recipient of the e-mail message,

$$N(t)a\mathrm{log}\left(\frac{t+c}{b}\right)n\left(t\right)=\frac{dN\left(t\right)}{dt}\frac{1}{t+c}.$$
(2)
Here the constant $`c`$ as a first approximation incorporates the fact that the measured time is not the true response time. The prime reason for the “shift” $`c`$ is that most people do not download new e-mail messages instantaneously but instead every 10 minutes or so. Furthermore, some time is obviously needed in order to formulate and write the reply. The values of the fit parameters are $`a=0.14`$, $`b=0.21`$ and $`c=0.25`$, see figure caption for details of the fitting procedure. We clearly see that equation (2) is an excellent approximation of the data over three decades. An additional feature of the data that can be rationalized is the first “ wiggle” occurring around $`1016`$ hours. Many people send e-mail messages just before leaving their work place. Since people generally share the same working hours (provided that they live in the same time zone), those messages are not answered before the next day. The main result to be extracted from this experiment compared to the previous one is that the response time distribution of the population exposed to new information is independent of the specific nature of the new information. This suggest that the observed behaviour indeed has its origin in psychological traits of the population. Furthermore, it should be noted that whereas the previous experiment deals with time scales of up to 2 months, this experiment deals with time scales up to only one week, i.e., a difference in time scales of approximately a factor of ten. Nevertheless, the same distribution is obtained. This suggest that the difference in time scales are mainly due to differences in communication channels, i.e., impersonal or personal.

## 4 Discussion

The results shown in figures 1 and 2 suggests that the general response rate of a human population $`n\left(t\right)`$ is that of a power law $`n\left(t\right)1/(t+c)`$ where $`c`$ depends on the details of the distribution channel(s) used. Such power law dependence of the rate of “events” are found in both natural and sociological systems. Two quite different examples are the Omori law for the rate of aftershocks as a function of time elapsed since the main shock  and the distribution of returns in the SP500 exceeding a specified threshold as a function of time elapsed since the “main event”, i.e., the “crash”, . A more appealing analogy, at least on a qualitative level, is provided by the relaxation of spin glasses subjected to a magnetic field: At time $`t`$ after the appearance of the interview (experiment I), the exposed population consists of two groups, namely those who have not downloaded a paper and those who have. Similarly with respect to experiment II, at any time $`t`$ the population considered consists of two groups, namely those who have an e-mail to answer and those who have not. The transition from the first state to the second demands the crossing of some threshold specific to each individual. We thus imagine that the announcement of the URL/the reception of emails plays the role a “field” to which the exposed population is subjected and study the relaxation process by monitoring the number of downloads/the number of replies as a function of time. Hence, we may view the process of downloading/replying as a diffusion process in a random potential, where the act of downloading/replying is similar to that of a barrier-crossing in the Trap model of spin glasses . The most pressing unanswered question raised by the two experiments is whether the observed power law only is a characteristic of the entire population or if it is also true for an individual over time. A qualitative argument that “ensemble averaging” may be the same as “time averaging” on the individual level (“ergodicity”) is that people react very differently to the same piece of information, i.e., one cannot generally deduce why people react on a specific piece of information. The reason is that individuals perceive the same piece of information quite differently due to individual psychology despite similar socio-economical status, e.g., not all rich people have an ambition to get richer because of the risk and/or work-load involved. This interpretations is supported by the results from the second experiment, where the new piece of information is changes with time and individual and nevertheless we get the same functional relationship between the time of the presentation of the new information and the response time distribution of the population. However, it seems a priori a quite formidable task to empirically verify whether these considerations are valid or not.

Acknowledgment: I thank D. Sergi, E. Moses and J.P. Eckmann for the data of experiment II.


text:
# Yet another galaxy identification for an ultraluminous X-ray sourceBased on observations collected at the Bologna Astronomical Observatory in Loiano, Italy, and at the Astronomical Observatory of Asiago, Italy.

## 1 Introduction

Ultraluminous X-ray sources (ULXs), i.e. off-nuclear sources with luminosities well above the Eddington limit for a typical neutron star, are now known to populate about 30% of nearby galaxies (see van der Marel 2003 for a review). Despite much effort, little is presently known about the nature of these sources, particularly since the counterparts at other wavelengths are difficult to find. Indeed, these sources are often superposed against regions of high surface brightness in their host galaxy. Thus, most of the optical identifications are obtained by spatial coincidence of the counterparts and comparison with broadband spectral characteristics of known stars in the host (Liu et al. 2002; Wu et al. 2002; Zezas et al. 2002; Zampieri et al. 2003). In other cases, it has been possible to study only the nearby environment in the ULX host (Pakull & Mirioni 2002; Wang 2002; Roberts et al. 2003). In one case (NGC4698-ULX1) the optical spectral features allowed a clear identification as a background BL Lac object (Foschini et al. 2002a).

Here we report the identification of the nature of another ULX. The counterparts in the infrared and optical bands were found, and spectroscopy revealed clear Balmer and forbidden transition emission lines. The derived redshift of $`z=0.217`$ indicates a background galaxy in this case also.

## 2 The source: NGC4168-ULX1

NGC 4168 is an E2 elliptical galaxy located in the Virgo cluster ($`d=16.8`$ Mpc). It hosts an active galactic nucleus (AGN), classified as a Seyfert 1.9 by Ho et al. (1997). The galaxy was observed on 4 December 2001 using the European Photon Imaging Camera (EPIC) on board the *XMM-Newton* satellite (see Foschini et al. 2002b; for the X-ray part, we refer in the following to the results obtained in this paper, unless explicitly stated). EPIC is composed of two instruments: the PN-CCD camera (Strüder et al. 2001) and two MOS-CCD detectors (Turner et al. 2001). The effective exposure time was $`17.4`$ ks.

One ULX was found apparently associated with NGC 4168, being inside the galaxy’s $`D_{25}`$ ellipse, at $`45^{\prime \prime }`$ from the optical centre of the galaxy. The ULX has coordinates (J2000) $`\alpha =12^\mathrm{h}12^\mathrm{m}14\stackrel{s}{.}5`$ and $`\delta =+13^{}12\mathrm{}48\mathrm{}`$, with an uncertainty radius of $`4^{\prime \prime }`$.

The counts were not sufficient to extract a spectrum, so we could only convert to physical units using the count rates derived from the *eboxdetect* task of the XMM-SAS software (v. 5.2). We found a count rate of $`4.0\pm 0.7`$ counts s-1 in the 0.5–10 keV band, which corresponds to a flux of $`(1.8\pm 0.3)\times 10^{14}`$ erg cm-2 s-1 adopting a conversion factor of $`3\times 10^{11}`$ counts cm2 erg-1, which in turn was derived by using a power-law model with $`\mathrm{\Gamma }=2.0`$ and an average Galactic column density of $`N_\mathrm{H}=3\times 10^{20}`$ cm-2. This value for $`\mathrm{\Gamma }`$ is common among the ULXs found with *XMM-Newton*. This flux value was corrected according to the energy encircled fraction (Ghizzardi 2001). Correction for vignetting (Lumb 2002) was not applied because the source is close to the center of the field of view ($`<`$$`2\mathrm{}`$). At the distance of NGC 4168 the resulting $`0.510`$ keV luminosity is $`6\times 10^{38}`$ erg s-1 (for a discussion on the luminosity threshold of ULXs, see Foschini et al. 2002c).

## 3 Archival optical and near-infrared data

A single optical object with coordinates (J2000) $`\alpha =12^\mathrm{h}12^\mathrm{m}14\stackrel{s}{.}60`$ and $`\delta =+13^{}12\mathrm{}47\stackrel{}{.}`$8 (i.e., 1$`\stackrel{}{.}`$4 from the *XMM-Newton* position, thus well within the X-ray error box) is clearly visible on the Digitized Sky Survey (DSS)111http://archive.eso.org/dss/dss/ observation of NGC 4168, originally made with the 48-inch Schmidt telescope at Palomar Observatory on 14 April 1955. We thus consider it as the optical counterpart of NGC4168-ULX1.

The source is also present in the US Naval Observatory (USNO) B1.0 Catalog (Monet et al. 2003) with the identification number $`10320222128`$. The magnitudes in the different bands are $`B1=19.1`$, $`B2=17.8`$, $`R1=18.1`$, $`R2=18.2`$, and $`I=17.8`$. The $`B1`$ and $`R1`$ magnitudes refer to the Palomar Observatory Sky Survey I (POSSI), performed between 1949 and 1965. The $`B2`$, $`R2`$ and $`I`$ magnitudes are measured from the Palomar Observatory Sky Survey II (POSSII), performed from 1985 to 2000.

Two Micron All Sky Survey (2MASS)222http://www.ipac.caltech.edu/2mass/ near-infrared (NIR) observations were performed on 5 April 2000, and the source was detected in all the three bands ($`JHK`$). From the Second Point Source Catalog, we get $`J`$ = 16.8$`\pm `$0.2, $`H`$ = 16.1$`\pm `$0.2, and $`K`$ = 15.3$`\pm `$0.1.

## 4 Optical observations at Loiano and Asiago

### 4.1 Spectra

Six medium-resolution optical spectra of the ULX were acquired between 21:23 UT of 4 March and 03:24 UT of 5 March 2003 in Loiano (Italy) with the Bologna Astronomical Observatory $`1.52`$ metre “G.D. Cassini” telescope plus BFOSC, for a total exposure time of 3 hours. One single 1-hour medium-resolution spectrum was also obtained at the “L. Rosino” Astronomical Observatory of Asiago (Italy), starting at 01:02 UT of 26 February 2003, with the $`1.82`$ metre “Copernicus” telescope plus AFOSC. The Cassini telescope was equipped with a $`1300\times 1340`$ pixels EEV CCD. The Copernicus telescope mounted instead a $`1100\times 1100`$ pixels SITe CCD. In both cases, Grism #4 and a slit width of $`2^{\prime \prime }`$ were used, providing nominal spectral coverages 3500-8500 Å (Loiano) and 4000-8000 Å (Asiago). The use of this setup secured a final dispersion of $`4.0`$ Å/pix and $`4.3`$ Å/pix for the spectra acquired in Loiano and Asiago, respectively.

Spectra, after correction for flat-field, bias and cosmic-ray rejection, were background subtracted and optimally extracted (Horne 1986) using IRAF333http://iraf.noao.edu/. Wavelength calibration of the Loiano spectra was performed using He-Ar lamps, while the one acquired in Asiago was calibrated with Cd-Hg-Ne lamps. Spectra taken in Loiano were then flux-calibrated by using the spectrophotometric standard Hiltner 600 (Hamuy et al. 1992, 1994) and finally stacked together. Correction for slit losses was also applied to the continuum. For the Asiago spectrum no spectroscopic standard was available. Wavelength calibration was checked by using the positions of background night sky lines; the error was $``$0.5 Å for both spectra.

### 4.2 Imaging

Between 23:45 UT of 4 March and 01:28 UT of 5 March 2003 $`UBVRI`$ photometry was also acquired in Loiano with BFOSC under an average seeing of 1$`\stackrel{}{.}`$5. The EEV CCD, with a scale of 0$`\stackrel{}{.}`$58/pix, secured a field of 12$`\stackrel{}{.}`$6$`\times `$12$`\stackrel{}{.}`$6. Images were corrected for bias and flat-field in the usual fashion and calibrated using the Rubin 149 field (Landolt 1992); the calibration accuracy is better than 3% in $`BVR`$ and better than 5% in $`UI`$. The source (in Fig. 1) is well detected in all bands.

Galactic absorption in the optical and NIR bands along the direction of the source was evaluated using the Galactic dust infrared maps by Schlegel et al. (1998); from these data we obtained a color excess $`E(BV)`$ = 0.036 mag. By applying the relation of Cardelli et al. (1989), we derived the following magnitudes for Galactic extinction correction: $`A_U`$ = 0.17, $`A_B`$ = 0.15, $`A_V`$ = 0.11, $`A_R`$ = 0.08, $`A_I`$ = 0.07, $`A_J`$ = 0.03, $`A_H`$ = 0.02, and $`A_K`$ = 0.01.

## 5 Results

The average spectrum taken at Loiano (Fig. 2) shows a number of emission features that can be readily identified with typical optical nebular lines. These include \[O ii\] $`\lambda `$3727, Hβ, \[O iii\] $`\lambda `$5007, Hα, \[N ii\] $`\lambda \lambda `$6548, 6583, and \[S ii\] $`\lambda \lambda `$6716, 6731. All identified emission lines yield a redshift of $`z=0.217`$. This result is, in hindsight, not surprising, as the images acquired in Loiano show that the field is well populated by galaxies of various sizes and brightnesses, with very few foreground stars. The Asiago spectrum, albeit of lower S/N due to its shorter exposure, shows the presence of Hβ and \[O ii\] $`\lambda `$3727 at the same redshift, thus confirming the result of the Loiano observations. Unfortunately, the actual useful spectral range secured by this spectrum has its red end at 7800 Å; therefore it does not cover the region containing Hα. Table 1 reports the emission-line fluxes as determined from the Loiano spectrum, dereddened for Galactic absorption. Given the limited S/N and resolution of the spectrum, no correction for starlight contamination (e.g., Ho et al. 1993, 1997) was attempted, but this does not strongly affect any of our conclusions.

Magnitudes of the source in the Loiano images were measured through aperture photometry, since its profile is significantly larger than the image PSF (2$`\stackrel{}{.}`$6 versus 1$`\stackrel{}{.}`$5) and, as measured on the $`U`$-band frame, is possibly elongated in the NW-SE direction. Using an aperture radius of 4.5 pixels (corresponding to 2$`\stackrel{}{.}`$6), we found the following optical magnitudes (not corrected for Galactic absorption): $`U`$ = 19.49$`\pm `$0.08, $`B`$ = 19.78$`\pm `$0.04, $`V`$ = 18.69$`\pm `$0.03, $`R`$ = 18.33$`\pm `$0.03, and $`I`$ = 17.68$`\pm `$0.05. The broad-band spectral energy distribution (SED) of the source, constructed with data from *XMM-Newton* (X-ray), Loiano (optical), and 2MASS (NIR), is shown in Fig. 3. The optical-NIR data points were corrected for Galactic absorption and converted into flux densities using the tables by Fukugita et al. (1995) for the optical and by Bersanelli et al. (1991) for the NIR. A 5% systematic error was added in quadrature to account for the uncertainties in the magnitude-to-flux conversion factors (Fukugita et al. 1995). No correction for any possible further absorption along the line of sight produced by the halo of NGC 4168 was included; this, however, should not be high as the source is well detected in the $`U`$ band.

## 6 Discussion

So, yet another putative ULX is found to be a background source. Assuming a cosmology with $`H_0`$ = 65 km s-1 Mpc-1, $`\mathrm{\Omega }_\mathrm{\Lambda }`$ = 0.7 and $`\mathrm{\Omega }_\mathrm{m}`$ = 0.3, we find that the luminosity distance to this source is $`d_L`$ = 1.16 Gpc, and that its X-ray luminosity is $`2.9\times 10^{42}`$ erg s-1 in the 0.6–12 keV rest-frame energy range. The angular size of the source translates into a linear diameter of about 30 kpc at $`z`$ = 0.217.

The measured value for the luminosity is $``$100–1000 times less that of “classical” active galactic nuclei (AGNs). Indeed, visual inspection of the optical spectrum in Fig. 2 suggests that the lines are due to stellar photoionization, rather than to an AGN. The diagnostic line ratios \[N ii\]/Hα, \[S ii\]/Hα, and \[O iii\]/Hβ, together with the nondetection of substantial \[O i\] $`\lambda `$6300 emission, place this source in the regime of metal-rich giant extragalactic H ii regions or starburst nuclei (Ho et al. 1993, 1997).

Inspection of the SED of the source (Fig. 3) shows that a peak in the optical-NIR domain is present, followed by a drop across the $`U`$ and $`B`$ bands. Given the dominance of the stellar continuum in the spectrum (Fig. 2), most of the optical and NIR light comes from the integrated stellar emission from the galaxy. Although the S/N of the spectrum is insufficient to place strong constraints on the age and metallicity of the stellar population, they are not inconsistent with that of an evolved population. The optical-NIR SED also supports this conjecture. The X-ray emission, by contrast, is likely to be associated with the star-forming regions that give rise to the optical line emission. This may consist of a single starburst nucleus or multiple off-nuclear sources. An X-ray luminosity of $``$1042 erg s-1 is high, but not unusual, for a starburst (David et al. 1992).

The strength of the optical emission lines, after accounting for internal reddening, can be used to estimate the star formation rate (SFR) and metallicity. Assuming an intrinsic Balmer decrement of Hα/Hβ = 2.86 (Osterbrock 1989) and the extinction law of Cardelli et al. (1989), the observed Hα/Hβ = 5.82 implies an internal reddening of $`E(BV)`$ = 0.72 mag. Following Kennicutt (1998), we determine a SFR of 43$`\pm `$5 $`M_{}`$ yr-1 from the reddening-corrected Hα luminosity of $`5.5\times 10^{42}`$ erg s-1. Similarly, the \[O ii\] luminosity yields a SFR of 80$`\pm `$20 $`M_{}`$ yr-1. Thus, to within a factor of $``$2 and within the uncertainties, the two estimates are consistent.

The optically derived SFR, on the other hand, is significantly lower than that predicted from the total X-ray luminosity. From the calibration of Ranalli et al. (2003), a 0.5–10 keV luminosity of 3$`\times `$1042 erg s-1, if entirely attributed to star formation, requires a SFR of $``$300 $`M_{}`$ yr-1, a factor of $``$5 higher than our estimates based on optical line emission. There are two possible explanations for this apparent discrepancy. First, the optically derived SFR may be an underestimate of the total SFR if star formation is more extended than the 2′′-wide slit: indeed, in case the distribution of star-forming regions is wider than that of the optical continuum emission, the correction for slit losses underestimates the actual total flux of emission lines. By contrast, the X-ray emission integrates over a larger area. Second, only a portion ($``$20%) of the X-ray emission originates from star formation, with the dominant fraction coming from an AGN. Although the optical spectrum shows no signs of nonstellar activity, we cannot exclude the presence of a low-level AGN, especially one that is heavily obscured. Both of these possibilities can be tested with further observations.

The detection of \[O ii\], \[O iii\], and Hβ also allows us to infer the gas-phase oxygen abundance. Following Kobulnicky et al. (1999), the $`R_{23}`$ parameter gives 12 + log (O/H) = 7.6 or 8.9. Considering the intrinsic luminosity of the source (rest-frame $`M_B`$ $``$21.5 mag) and its \[O iii\]/\[N ii\] ratio, the larger of the two values is likely to be the correct one (see Kobulnicky et al. 1999), pointing to a basically solar oxygen abundance.

The comparison of the USNO archival magnitudes with the results of the optical photometry acquired in Loiano seems to suggest long-term variability for the source. However, we caution the reader that it is not infrequent to find discrepancies as large as 0.6 mag (see e.g. Masetti et al. 2003) between the automatically extracted archival magnitudes and the actual ones, calibrated through Landolt (1992) fields. This is most likely due to the fact that magnitude extraction from archival images was performed by using automatic pipelines and was (for the oldest surveys) applied to photographic plates.

In the minisurvey in Foschini et al. (2002b) 18 ULXs in 10 nearby galaxies were detected. According to the statistics from the Lockman hole studies (Hasinger et al. 2001) we expect to find 3.2 background sources with 2–10 keV flux higher than $`10^{14}`$ erg cm-2 s-1 in an area corresponding to the sum of all the $`D_{25}`$ ellipses of the considered galaxies. If we also take into account the case reported in Foschini et al. (2002a), with the present identification we already found 2 background objects among the 18 ULX candidates: so, we are close to the statistical limit by assuming that the Lockman Hole studies are valid over a wide range of the sky. Albeit the two selected cases might suffer from an observational bias (they are quite bright and outside the main body of the parent galaxy), additional detection of other background objects in the minisurvey should be considered with care.

In conclusion, we wish to remark that this case, as that of the ULX in NGC 4698 reported by Foschini et al. (2002a), has been instructive: despite the availability of broad-band multiwavelength data, they have proven insufficient to disclose the true nature of the source. The most crucial piece of information is optical spectroscopy of sufficient quality to reveal clear spectral features, which simultaneously provide redshift information and spectral diagnostics.

###### Acknowledgements.

We wish to thank L. Zampieri for the coordination of the Asiago Service Mode observations, and P. Grandi and M. Cappi for useful discussions. We also thank R. Gualandi for the night assistance in Loiano and H. Navasardyan for having performed the observations in Asiago as a Service Mode run. This research has made use of the NASA Astrophysics Data System Abstract Service, of the NASA/IPAC Extragalactic Database (NED), and of the NASA/IPAC Infrared Science Archive, which are operated by the Jet Propulsion Laboratory, California Institute of Technology, under contract with the National Aeronautics and Space Administration.


text:
# Minimal model for beta relaxation in viscous liquids

## Abstract

Contrasts between beta relaxation in equilibrium viscous liquids and glasses are rationalized in terms of a double-well potential model with structure-dependent asymmetry, assuming structure is described by a single order parameter. The model is tested for tripropylene glycol where it accounts for the hysteresis of the dielectric beta loss peak frequency and magnitude during cooling and reheating through the glass transition.

Viscous liquids approaching the calorimetric glass transition have extremely long relaxation times viscliq . The main relaxation is termed the alpha relaxation. There is usually an additional minor “beta” process at higher frequencies. Dielectric relaxation is a standard method for probing liquid dynamics diel . The study of dielectric beta relaxation in simple viscous liquids was pioneered by Johari and Goldstein more than 30 years ago joh70 ; joh85 , but the origin of beta relaxation is still disputed kah97 ; wag98 ; joh02a . It is unknown whether every molecule contributes to the relaxation all  or only those within “islands of mobility” gol69 ; joh76 ; kop00 . Similarly, it is not known whether small angle jumps all ; kau90 ; vog00  or large angle jumps arb96  are responsible for the beta process.

Improvements of experimental techniques have recently lead to several new findings. The suggestion ols98 ; leo99 ; wag99  that the excess wing of the alpha relaxation usually found at high frequencies is due to an underlying low-frequency beta process was confirmed by long time annealing experiments by Lunkenheimer and coworkers lunken  (an alternative view is that the wing is a non-beta type relaxation process paluch ). This lead to a simple picture of the alpha process: Once the effect of interfering beta relaxation is eliminated, alpha relaxation obeys time-temperature superposition with a high frequency loss $`\omega ^{1/2}`$ ols01 . Moreover, it now appears likely that all liquids have one or more beta relaxations lunken ; han97 ; kud97 ; joh02b : Liquids like propylene carbonate, glycerol, salol, and toluene are now known to possess beta relaxation, while o-terphenyl, previously thought to have a beta process only in the glassy state, has one in the equilibrium liquid phase as well. Finally, it has been shown that beta relaxation in the equilibrium liquid does not behave as expected by extrapolation from the glassy phase: In some cases the beta loss peak frequency is temperature independent in the liquid phase (e.g., sorbitol ols98 ), in other cases it is very weakly temperature dependent. On the other hand, the beta relaxation strength always increases strongly with temperature in the liquid phase ols00 ; note1 .

The contrasts between beta relaxation in liquid and glass are clear from Fig. 1(a) which shows beta loss peak frequency and maximum loss for tripropylene glycol cooled through the glass transition and subsequently reheated note2 .

In the glassy phase (at low temperatures) the loss peak frequency is strongly temperature dependent while the maximum loss varies little. On the other hand, the loss is strongly temperature dependent in the equilibrium liquid phase. Here we even see the loss peak frequency decreasing upon heating. How is one to understand these findings? A clue is provided by Fig. 1(b) which shows the maximum loss as a function of time after an “instantaneous” temperature step, i.e., instantaneous on the time scale of structural (alpha) relaxation. This experiment utilizes a special-purpose setup with a cell consisting of two aluminum discs separated by three capton spacers (layer distance $`20\mu `$, empty capacitance 30 pF). One disc, where temperature is measured via an NTC resistor, is placed on a Peltier element. Less than six seconds after a 0.672 K temperature jump is initiated, temperature is stable within 1 mK. In this setup we measure at 10 kHz which is the loss peak frequency (changes of loss peak frequency lead only to second order corrections of $`ϵ_{\mathrm{max}}^{\prime \prime }`$). The sampling time is 2 s. Figure 1(b) shows a very fast change of the maximum loss, followed by relaxation toward the equilibrium value taking place on the structural (alpha) relaxation time scale. The existence of an instantaneous increase of the loss clearly indicates a pronounced asymmetry of the relaxing entity. Inspired by this fact we adopt the standard asymmetric double-well potential model (Fig. 2) with transitions between the two free energy minima described by rate theory diel ; glref ; pol72 ; gil81 ; buc01 .

In terms of the small barrier $`U`$ and the asymmetry $`\mathrm{\Delta }`$, loss peak frequency $`f_{\mathrm{max}}`$ and maximum loss $`ϵ_{\mathrm{max}}^{\prime \prime }`$ are given gil81  by

$$\{\begin{array}{ccc}f_{\mathrm{max}}& =& f_0\mathrm{exp}\left(\frac{2U+\mathrm{\Delta }}{2k_BT}\right)\mathrm{cosh}\left(\frac{\mathrm{\Delta }}{2k_BT}\right)\hfill \\ ϵ_{\mathrm{max}}^{\prime \prime }& =& ϵ_0^{\prime \prime }(T)\mathrm{cosh}^2\left(\frac{\mathrm{\Delta }}{2k_BT}\right)\hfill \end{array}.$$
(1)
The prefactor $`f_0`$ is assumed to be structure and temperature independent while $`ϵ_0^{\prime \prime }(T)1/T`$ gil81  is assumed to be structure independent: $`ϵ_0^{\prime \prime }(T)=T_0/T`$. The free energy differences $`U`$ and $`\mathrm{\Delta }`$ are expected to change with changing structure, but freeze at the glass transition. In terms of the fictive temperature $`T_f`$ our model is based on

$$\{\begin{array}{ccccc}U& =& U_0& +& ak_BT_f\\ \mathrm{\Delta }& =& \mathrm{\Delta }_0& & bk_BT_f\end{array}.$$
(2)
Equation (2) follows from minimal assumptions: Suppose structure is parameterized by just one variable, $`s`$. Only a rather narrow range of temperatures is involved in studies of beta relaxation in the liquid phase and around the glass transition. Consequently, structure varies only little and, e.g., $`U(s)`$ may be expanded to first order: $`U(s)=c_0+c_1(ss_0)`$. For the equilibrium liquid, $`s=s(T)`$ which may also be expanded to first order. By redefining $`s`$ via a linear transformation we obtain $`s=T`$ at equilibrium while $`U(s)`$ is still linear in $`s`$. A single variable describing structure, which at equilibrium is equal to temperature, is – consistent with Tool’s 1946 definition too46  – to be identified with the fictive temperature: $`s=T_f`$. Thus one is lead to Eq. (2) where signs are chosen simply to ensure $`a,b>0`$ in fit to data.

The model has 6 parameters: $`f_0`$, $`T_0`$, $`U_0`$, $`\mathrm{\Delta }_0`$, $`a`$, and $`b`$. These were determined by measuring the instantaneous changes of $`f_{\mathrm{max}}`$ and $`ϵ_{\mathrm{max}}^{\prime \prime }`$ upon a temperature step, as well as their thermal equilibrium changes. Figure 1(b) allows determination of the instantaneous change of the loss. It is not possible to determine the instantaneous change of $`f_{\mathrm{max}}`$. Instead we extrapolated measurements obtained by the standard cell note2  as follows (Fig. 3): Beta loss is monitored by first annealing at 183.0 K, subsequently changing temperature to 181.0 K. The latter data show a linear relation between $`\mathrm{log}f_{\mathrm{max}}`$ and $`\mathrm{log}ϵ_{\mathrm{max}}^{\prime \prime }`$ which, knowing the instantaneous change per Kelvin of $`ϵ_{\mathrm{max}}^{\prime \prime }`$ from Fig. 1(b), is extrapolated to short times.

In the data analysis it is convenient to eliminate $`\mathrm{cosh}`$ by introducing the variable

$$Y\left(\frac{f_{\mathrm{max}}}{f_0}\right)^2\left(\frac{ϵ_{\mathrm{max}}^{\prime \prime }}{ϵ_0^{\prime \prime }(T)}\right).$$
(3)
On the fast time scale $`U`$ and $`\mathrm{\Delta }`$ are frozen so Eq. (1) implies

$$\{\begin{array}{ccc}\frac{d\mathrm{ln}Y}{d\mathrm{ln}T}|_{\mathrm{inst}}& =& \frac{2U+\mathrm{\Delta }}{k_BT}\\ \frac{d\mathrm{ln}ϵ_{\mathrm{max}}^{\prime \prime }}{d\mathrm{ln}T}|_{\mathrm{inst}}& =& \frac{\mathrm{\Delta }}{k_BT}\sqrt{1\frac{ϵ_{\mathrm{max}}^{\prime \prime }}{ϵ_0^{\prime \prime }(T)}}1\end{array}.$$
(4)
Having determined the instantaneous changes of $`ϵ_{\mathrm{max}}^{\prime \prime }`$ and $`f_{\mathrm{max}}`$, Eqs. (1) and (4) provide 4 equations for the 6 model parameters. The last two equations come from the temperature dependence of loss and loss peak frequency at thermal equilibrium where $`U(T)`$ and $`\mathrm{\Delta }(T)`$ are given by Eq. (2) with $`T_f=T`$, leading to

$$\{\begin{array}{ccc}\frac{d\mathrm{ln}Y}{d\mathrm{ln}T}|_{\mathrm{eq}}& =& \frac{2U_0+\mathrm{\Delta }_0}{k_BT}\\ \frac{d\mathrm{ln}ϵ_{\mathrm{max}}^{\prime \prime }}{d\mathrm{ln}T}|_{\mathrm{eq}}& =& \frac{\mathrm{\Delta }_0}{k_BT}\sqrt{1\frac{ϵ_{\mathrm{max}}^{\prime \prime }}{ϵ_0^{\prime \prime }(T)}}1\end{array}.$$
(5)
Using Eq. (4) for the data of Figs. 1(b) and 3, and using Eqs. (1) and (5) for equilibrium state measurements at 183 and 185 K, the 6 parameters determined for tripropylene glycol note3  are: $`f_0=6.210^{11}`$ Hz, $`T_0=130`$ K, $`U_0/k_B=1107`$ K, $`\mathrm{\Delta }_0/k_B=3039`$ K, $`a=23.28`$, $`b=12.77`$. $`U_0<0`$ reflects the fact that the beta loss peak frequency in the liquid phase decreases as temperature increases \[Fig. 1(a)\]. Physically, this anomalous behavior is caused by the barrier $`U`$ increasing more than $`T`$ upon heating.

Once all parameters are fixed the model predicts how $`f_{\mathrm{max}}`$ and $`ϵ_{\mathrm{max}}^{\prime \prime }`$ correlate for the continuous passage through the glass transition of Fig. 1(a). To analyze these data within the model we first note that if $`y(x)=1/\mathrm{cosh}^2(x)`$ then $`x=\varphi (y)`$ where $`\varphi (y)=\mathrm{ln}(1/\sqrt{y}+\sqrt{1/y1})`$, so Eq. (1) may be inverted: $`\mathrm{\Delta }/2k_BT=\varphi (ϵ_{\mathrm{max}}^{\prime \prime }/ϵ_0^{\prime \prime }(T))\mathrm{\Phi }`$. Since $`\mathrm{ln}Y=(2U+\mathrm{\Delta })/k_BT`$ both $`\mathrm{\Phi }`$ and $`\mathrm{ln}Y`$ involve fictive temperature, and consequently both exhibit hysteresis at the glass transition. However, fictive temperature is eliminated by considering the following variable

$`Z(f_{\mathrm{max}},ϵ_{\mathrm{max}}^{\prime \prime })`$ $``$ $`\left(2ab\right)\mathrm{\Phi }{\displaystyle \frac{b}{2}}\mathrm{ln}Y`$ (6)
$`=`$ $`{\displaystyle \frac{a\mathrm{\Delta }_0+bU_0}{k_BT}}.`$
The model is tested in Fig. 4. The glass transition is not visible and most hysteresis is eliminated (better elimination was obtained in Ref. ols00  but without theoretical basis and with one free parameter note4 ). The line shown is the prediction of Eq. (6).

In conclusion, a minimal model for beta relaxation in viscous liquids has been proposed. The model is built on the four simplest possible assumptions: 1) Beta relaxation involves only two levels, 2) structure is determined by just one order-parameter, 3) first order Taylor expansions apply, 4) the two characteristic free energies $`U`$ and $`\mathrm{\Delta }`$ freeze at the glass transition. The model is clearly oversimplified. For instance, it predicts a Debye response which is not observed, and $`U`$ and $`\mathrm{\Delta }`$ would be expected to vary slightly with temperature in the glass. Nevertheless, the model is able to rationalize the contrasts between beta relaxation in liquids and in glasses. One final puzzling observation should be mentioned: The asymmetry $`\mathrm{\Delta }`$ extrapolates to zero at a temperature which is close to the temperature where alpha and beta relaxations merge. We have seen the same phenomenon in sorbitol, a pyridine-toluene solution, polypropylene-glycol-425, and in 4,7,10-trioxydecane-1,13-diamine ols98 ; ols00 , and have found no exceptions. This finding indicates that the merging temperature is fundamental, a symmetry is somehow broken below this temperature.

###### Acknowledgements.

This work was supported by the Danish Natural Science Research Council.