text:
# Average lattice symmetry and nanoscale structural correlations in magnetoresistive manganites

\[

## Abstract

We report x-ray scattering studies of nanoscale structural correlations in the paramagnetic phases of the perovskite manganites La0.75(Ca0.45Sr0.55)0.25MnO3, La0.625Sr0.375MnO3, and Nd0.45Sr0.55MnO3. We find that these correlations are present in the orthorhombic $`O`$ phase in La0.75(Ca0.45Sr0.55)0.25MnO3, but they disappear abruptly at the orthorhombic-to-rhombohedral transition in this compound. The orthorhombic phase exhibits increased electrical resistivity and reduced ferromagnetic coupling, in agreement with the association of the nanoscale correlations with insulating regions. In contrast, the correlations were not detected in the two other compounds, which exhibit rhombohedral and tetragonal phases. Based on these results, as well as on previously published work, we propose that the local structure of the paramagnetic phase correlates strongly with the average lattice symmetry, and that the nanoscale correlations are an important factor distinguishing the insulating and the metallic phases in these compounds.

. \]

The physical mechanism underlying the magnetic-field-induced insulator-metal transition in perovskite manganites $`A_{1x}B_x`$MnO3 has been the subject of intense experimental and theoretical investigation since its rediscovery in 1993 . One of the motivations for such considerable attention is the unusually large diminution of the electrical resistivity observed at the magnetic-field-induced transition, which is now commonly referred to as the colossal magnetoresistance (CMR). Several different kinds of the CMR effect are known . The most widely studied variant of this effect is the transition from a paramagnetic insulating (PI) to a ferromagnetic metallic (FM) phase. The large difference between the resistivities of these two phases lies at the heart of the CMR effect. The metallic nature of the FM phase has been explained within the framework of the double-exchange mechanism . However, the physical mechanism responsible for the large resistivity of the PI phase remains poorly understood.

The situation is complicated significantly by the fact that transport properties of the paramagnetic state in manganites with the same doping level $`x`$, but with different cations $`A`$ and $`B`$, are often very different. For example, the electrical resistivity of La0.7Sr0.3MnO3 shows a metallic behavior (i.e. grows with temperature) , while the resistivity of La0.7Ca0.3MnO3 exhibits the temperature dependence typical of an insulator, decreasing with temperature . As a consequence, the resistivity of the latter compound is significantly larger than that of the former. These differences cannot be explained by steric modification of the electronic bandwidth $`W`$ due to the evolution of the average structural parameters from one composition to another as a result of the differing size of the divalent dopants . In fact, the value of $`W`$ in La0.7Ca0.3MnO3 is expected to differ from that in La0.7Sr0.3MnO3 by less than 1% as a result of such effects . The widely accepted solution to this problem, first proposed theoretically , is that the electron-phonon coupling $`\lambda `$ plays an essential role. The resistivity is expected to be larger in compounds with large values of $`\lambda `$.

There is a significant amount of experimental evidence that electron-phonon coupling plays a key role in the manganites. It is well established that small lattice polarons are present in the PI state . A lattice polaron forms when an $`e_g`$ electron localizes on a Mn3+ ion, and the surrounding oxygen octahedron distorts due to the Jahn-Teller effect. Formation of lattice polarons leads to an increase of the electrical resistivity. In addition to single polarons, the paramagnetic state exhibits correlated lattice distortions with a correlation length of several lattice constants . It has been proposed that the structural correlations reported in Ref.  reflect the presence of nanoscale regions possessing charge and orbital order, and that these regions are responsible for the enhanced resistivity of the PI phase. The temperature dependence of the electrical resistivity has been shown to follow the concentration of the correlated regions , in agreement with this hypothesis. Further, the correlations are strongly suppressed in an applied magnetic field , implying that they play an important role in the CMR effect. The detailed structure of the correlated regions remains to be elucidated , but recent experimental studies suggest  that the correlated regions exhibit CE-type and striped charge and orbital ordered structures similar to those found in manganites possessing long-range charge and orbital order . Numerous theoretical works have also been devoted to the problem of structural and electronic inhomogeneities in manganites , but no realistic description of the PI state has been achieved thus far.

In the above scenario, the differences in the transport properties of the PI state in various manganites must then be explained by variations in the concentration, size, or other properties of the correlated regions. Such variations might, in principle, result from a dependence of the electron-phonon coupling $`\lambda `$ on sample composition. However, $`\lambda `$ is not expected to depend significantly on the $`A`$-site ionic composition (because it is, essentially, a local parameter proportional to the Jahn-Teller energy $`E_{JT}`$). Thus, other reasons for the observed variations of physical properties must be identified. Among the possible explanations, it has been proposed, for example, that in manganites with a large variance in the $`A`$ cation radius, cation size disorder could produce strain fields stabilizing local Jahn-Teller distortions . This explanation does not work, however, in the case of the La0.7Sr0.3MnO3 and La0.7Ca0.3MnO3 compounds discussed above, since the metallic Sr compound exhibits larger A-site cation radius variation than the Ca compound, which is insulating at high temperatures.

In this paper, we discuss a different scenario. Specifically, we point out that the symmetry of the average structure can affect the energetics for the formation of local distortions by allowing certain types of distortions to occur more easily, or, alternatively, by suppressing them. The importance of the average long-range structure stems from the fact that the significant lattice distortions associated with the correlated regions give rise to long-range anisotropic strain . Thus, lattice symmetry is one of the important factors determining how easily the local distortions are accommodated in a given crystal structure. Manganites exhibit structural transitions between distorted perovskite phases possessing orthorhombic, rhombohedral, tetragonal, and monoclinic lattice symmetry. As explained above, it is reasonable to expect that some of these phases will be more susceptible towards formation of local distortions than others.

In this work, therefore, we investigate the connection between the average lattice symmetry and the local (nanoscale) correlations in several perovskite manganite compounds. X-ray scattering studies of the structural correlations in manganites exhibiting orthorhombic, rhombohedral, and tetragonal structural phases are reported. We find that in the samples investigated, as well as in the previously studied manganites, these correlations are present only in the orthorhombic $`O`$ phase. This result is most convincingly demonstrated in experiments with La0.75(Ca0.45Sr0.55)0.25MnO3. This compound undergoes an orthorhombic-to-rhombohedral transition with increasing temperature. We find that the structural correlations abruptly disappear at the transition to the rhombohedral state. Based on these results and on previously published work, we propose that the correlations form most easily in the orthorhombic $`O`$ phase. In this phase, the crystal lattice is contracted along the $`c`$ axis, as is also the case in structures which exhibit charge-ordered phases, in which the long axes of the Mn3+O6 octahedra lie in the $`ab`$ plane. In contrast, the correlations are suppressed in rhombohedral and tetragonal phases, in which the MnO6 octahedra are either undistorted, or elongated along the $`c`$ axis, on average. Thus, the local structure of the paramagnetic phase and, consequently, its physical properties do, in fact, correlate strongly with the average lattice symmetry.

Single crystals of La0.75(Ca0.45Sr0.55)0.25MnO3 (LCSMO), La0.625Sr0.375MnO3, and Nd0.45Sr0.55MnO3 were grown using the floating zone technique. X-ray diffraction measurements were carried out at beamline X22C at the National Synchrotron Light Source. A 10.3-keV x-ray beam was focused by a mirror, monochromatized by a double-crystal Ge (111) monochromator, scattered from the sample, and analyzed with a pyrolytic graphite crystal. The samples were mounted in a closed-cycle refrigerator (T=6-450 K). In this paper, Bragg peaks are indexed in the orthorhombic $`Pbnm`$ notation, in which the $`a`$, $`b`$, and $`c`$ axes run along the (1,1,0), (1,-1,0), and (0,0,1) cubic perovskite directions, respectively. Scattering vectors ($`h,k,l`$) are given in reciprocal lattice units.

We first focus our discussion on the properties of the LCSMO sample. Its high temperature PI state exhibits rhombohedral symmetry. With decreasing temperature, this compound first undergoes a first order structural transition into the orthorhombic PI phase at $`T_s`$360 K, and then becomes a ferromagnetic metal at $`T_c`$300 K. These transitions are clearly reflected in the behavior of the electrical resistivity, shown as solid line in Fig. 1. The x-ray scattering patterns around the (4,4,0) Bragg peak in the orthorhombic ($`T`$=300 K) and rhombohedral ($`T`$=400 K) phases are shown in Fig. 2. The sharp peaks at (4.5,3.5,0) and (3.5,4.5,0) at $`T`$=300 K arise from crystallographic twinning in the orthorhombic phase. The “butterfly”-shaped feature in the center is commonly attributed to scattering from uncorrelated polarons, also known as Huang scattering, together with thermal-diffuse scattering . At $`T`$=300 K, there are also 4 broad peaks located at (4$`\pm `$0.5,4,0) and (4,4$`\pm `$0.5,0). These peaks are best seen in scans along the $`h`$ or $`k`$ directions. One such scan at $`T`$=320 K is shown in Fig. 3(a).

The broad peaks at (4$`\pm `$0.5,4,0) and (4,4$`\pm `$0.5,0) reflect the presence of nanoscale correlated regions which are believed to possess charge and orbital order . The intensity of these peaks reflects the concentration of the correlated regions . To study the temperature dependence of this concentration, we took data at different temperatures. Some of these data are shown in Fig. 3(a). The data were fitted to a sum of a Gaussian line shape and a monotonically sloping background, the latter described by a power-law function. This background contains contributions from single polaron and thermal-diffuse scattering . The error bars were estimated as described in Ref. .

The intensity of the correlated peak in LCSMO as a function of temperature is shown in Fig. 4(a). The central observation of this work is that the correlated regions abruptly disappear at the transition from the orthorhombic to rhombohedral phase. Thus, in LCSMO, the structural correlations are present only in the orthorhombic state. The data of Fig. 4(a) show that the correlated regions also disappear in the FM state, as was reported in earlier works . However, in contrast to the FM state, the rhombohedral phase is insulating. Therefore, the suppression of the correlated insulating regions in this phase is a much more surprising result.

Fig. 4(b) shows that the electrical resistivity abruptly decreases and the magnetic susceptibility increases as the sample enters the rhombohedral phase and the correlated regions disappear. The Curie temperature and the effective magnetic moment extracted from the susceptibility data are $`T_c`$=302(2) K, $`p`$=6.2, and $`T_c`$=322(3) K, $`p`$=5.6 for the orthorhombic and the rhombohedral states, respectively. Thus, the orthorhombic state possesses the larger effective moment, but smaller ferromagnetic exchange constant than the corresponding values in the rhombohedral state. This observation is in agreement with the proposed CE-type charge and orbital ordered structure of the correlated regions, which should be insulating and exhibit antiferromagnetic correlations, provided that these regions are sufficiently large . We note, however, that neutron scattering experiments have thus far failed to detect the CE-type antiferromagnetic correlations in the PI state , and therefore the character of the magnetic interaction in the correlated regions is yet to be established experimentally. Alternatively, the increased effective ferromagnetic exchange constant of the more conducting rhombohedral phase can stem from the effects of double exchange which favors ferromagnetic correlations.

The ratio of the integrated intensity of the (4,4.5,0) peak to the integrated intensity of the main Bragg (4,4,0) peak in LCSMO is of the order of 10-5. In the long-range charge-ordered phases, such as that in La0.5Ca0.5MnO3, the typical intensity ratio between the brightest charge-ordering peaks and Bragg peaks is of the order of 10-2 (Ref. ). This is significantly larger than what is observed in LCSMO. One can estimate the intensity ratio of the (4,4.5,0) and the (4,4,0) peaks using the structural models proposed in Ref.  (the CE-type charge and orbital ordered structure) for La0.5Ca0.5MnO3, and in Ref.  (the ordering of Zener polarons) for Pr0.6Ca0.4MnO3. The calculated ratios are 2.7$`\times `$10-2 for the former, and 4.1$`\times `$10-2 for the latter model structure. Assuming that the correlated regions in LCSMO are small charge ordered regions with the structure identical to the charge-ordered state in La0.5Ca0.5MnO3, we obtain the result that the volume fraction of these regions (defined here as the fraction of MnO6 octahedra experiencing the correlated distortion) is smaller than 0.1%. In contrast, a much larger estimate of “a few percent” was obtained for the same fraction in a recent electron diffraction study of La2/3Ca1/3MnO3 . The reason for such a discrepancy is currently unclear. In either case, it is hard to see how such an apparently small volume fraction of the charge and orbital ordered insulating phase could significantly affect the transport properties.

Thus, these observations appear to contradict the strong correlation between the nanoscale correlations and the transport properties observed in a large number of manganites of different chemical compositions and doping levels . One way to explain this contradiction is to assume that the lattice distortions are significantly (at least a factor of 10) smaller in the correlated regions than they are in the corresponding long-range ordered phases. Since the diffraction signal from the correlated regions is proportional to the square of the magnitude of the lattice distortion, much larger volume fractions of the charge and orbital ordered insulating phase are then deduced from the experimental data under this assumption.

Alternatively, more elaborate models of the inhomogeneous PI state should be proposed. In particular, dynamic inhomogeneities might play an important role in this state, as evidenced by recent experiments in the related bi-layered manganites , in which the structural correlations were found to become completely dynamic at high temperatures. Also, the role of the uncorrelated local distortions must be clarified. The latter distortions are known to be present in both the rhombohedral and the orthorhombic manganites, and to have composition-independent distortion magnitude . In the orthorhombic La1-xCaxMnO3 system, the number of MnO6 octahedra experiencing a distortion (correlated or otherwise) was found to be functionally correlated to the magnetic and transport properties . The relative contributions of the correlated and uncorrelated distortions to these functional relationships remain to be characterized. It also remains to be seen whether these relationships hold on the rhombohedral side of the structural transition boundary. In this work, we do not attempt to characterize the uncorrelated distortions because it is difficult to distinguish scattering due to these distortions from thermal-diffuse scattering in our energy-integrating measurements. Thus, while it is clear that the correlations play an important role in the manganites, more work is needed to elucidate the microscopic structure of the PI state.

We next briefly discuss the other samples. La0.625Sr0.375MnO3 is rhombohedral (space group $`R\overline{3}c`$) at all temperatures and exhibits a paramagnetic state for $`T>`$370 K. The paramagnetic insulating state of Nd0.45Sr0.55MnO3 is tetragonal (space group $`I4/mcm`$) . In the latter sample, the PI state is observed for $`T>T_N`$220 K. The electrical resistivities of these samples are shown in Fig. 1. In the both samples, we did not detect any evidence of the structural correlations at any temperature. Figs. 3(b,c) show some of the collected data for these samples (Bragg peaks are again indexed using the orthorhombic notation). In contrast, these correlations are clearly present in the paramagnetic orthorhombic phase of La0.75(Ca0.45Sr0.55)0.25MnO3 discussed above, and also in the PI phases of numerous other orthorhombic manganites, such as Nd1-xSrxMnO3 (0.3$`x`$0.5) , La1-xCaxMnO3 (0.2$`x`$0.3) , and Pr1-xCaxMnO3 (0.3$`x`$0.5) . We suggest then that the nanoscale structural correlations described above can only occur in an orthorhombic phase. Of course, more systematic studies are needed to verify this intriguing hypothesis, but even absent such studies, we feel that there is sufficient evidence in hand to warrant speculation as to why such an observation might be true.

All the orthorhombic samples discussed above exhibit the well known $`O`$ structure with the lattice constants $`c/\sqrt{2}<ba`$ . An orthorhombic state with $`c/\sqrt{2}<ba`$ is also realized in the manganites with long-range CE-type charge and orbital order. In this structure, the $`c`$-axis contraction reflects the global effect of ordered Jahn-Teller distortions of the Mn3+O6 octahedra with the long axis lying in the $`ab`$ plane . Thus, on average, the long axes of the MnO6 octahedra lie in the $`ab`$ plane in both the long-range charge ordered phase, and in the $`O`$ paramagnetic insulating state . In contrast, the MnO6 octahedra are undistorted in the rhombohedral $`R\overline{3}c`$ state , and are significantly elongated along the $`c`$-axis in the tetragonal $`I4/mcm`$ state .

It is not completely surprising that the MnO6 octahedra are distorted in a similar manner in the manganites exhibiting nanoscale charge-ordered regions and in the manganites with long-range charge and orbital order. Note, however, that this is a quite nontrivial result because there is no a priori requirement that the local symmetry of the nanoscale ordered regions should match the average symmetry of the crystal lattice. The energetically favorable response of the orthorhombic structure to the long-range strains produced by the local distortions is likely to be one of the key factors explaining the reported observations. Theoretical calculations taking into account the details of the lattice structure are needed to confirm these ideas. In any case, the proposed correlation between the presence of the local lattice distortions and the average crystallographic structure is intriguing and deserves further theoretical and experimental investigation.

Finally, we return to the problem of the drastic variation of the electrical properties of the paramagnetic phase with chemical composition. We emphasize that here we are discussing the high-temperature paramagnetic phase, and not the low-temperature state. Specifically, we will consider the La1-x(Ca1-ySry)xMnO3 series of compounds with $`x`$ in the vicinity of 0.3. This system exhibits a well defined insulating behavior for $`y`$=0 (see Fig. 1), and is metallic for $`y`$=1 (Fig. 1, and Ref. ). At fixed $`x`$ and with $`y`$ increasing from 0 to 1, the system should, therefore, undergo an insulator-to-metal transition. It is well known that with $`y`$ increasing from 0 to 1, this system also undergoes an orthorhombic-to-rhombohedral structural transition . With the exception of the narrow region in the vicinity of the structural transition, the electrical resistivity shows the temperature dependence typical for an insulator in the orthorhombic samples, and that of a metal in the rhombohedral state . These data show that the metal-insulator transition coincides approximately with the structural transition in this series of compounds.

In the vicinity of the structural transition, the electrical resistivity exhibits temperature-dependent behavior typical of an insulator in the both structural phases, see Fig. 4(b), and Ref. . We would like to point out an interesting possibility that the insulating behavior of the rhombohedral phase could result from the effects of chemical and magnetic disorder. These effects are, in fact, well known to play an important role in mixed-valence manganites . In the La1-x(Ca1-ySry)xMnO3 series of compounds, in particular, the importance of disorder is evidenced by the large magnitude of the electrical resistivity deep in the paramagnetic metallic region of the phase diagram . It is, therefore, reasonable to expect that in samples with less chemical disorder, the metal-insulator transition boundary lies closer or even coincides with the structural transition boundary. We, therefore, believe that further theoretical and experimental investigation of the role of disorder in the La1-x(Ca1-ySry)xMnO3 series of compounds is of great interest.

The experimental data for one of the compounds from this series, La0.75(Ca0.45Sr0.55)0.25MnO3, were discussed in detail above. Our experiments clearly show that the nanoscale structural correlations are present only on the orthorhombic side of the phase boundary. In addition, these correlations are clearly present in the $`y`$=0 compound La0.75Ca0.25MnO3 , and are undetectable in the $`y`$=1 material La0.625Sr0.375MnO3. Thus, the existing experimental data show that the metallic rhombohedral state in this compound series does not exhibit the nanoscale correlations, while the correlations are present in the insulating orthorhombic state. The presence of the nanoscale insulating regions, therefore, appears to be an important factor distinguishing the insulating and the metallic phases. Thus, the nanoscale correlations should be included in any physical mechanism explaining the observed large variation of the high-temperature properties of manganites and, consequently, the observed values of the CMR effect.

In summary, we report x-ray scattering studies of nanoscale structural correlations in several manganite compounds exhibiting orthorhombic, rhombohedral, and tetragonal perovskite phases. Based on these results, and also on previously published work, we propose that these correlations form most easily in the orthorhombic $`O`$ phase, and are suppressed in the rhombohedral and the tetragonal phases. Thus, the local structure of the paramagnetic phase appears to correlate strongly with the average lattice symmetry. Our experiments strongly suggest that changes in physical properties observed at structural transitions in manganites cannot be understood through spatially uniform phases, and that nanoscale inhomogeneities must play an essential role in these transitions.

Note added: A recent work by Mira at al. reports measurements of thermal properties, magnetoresistance, and lattice expansion in various manganites . The analysis of these quantities led the authors of Ref.  to a conclusion that lattice effects (electron-lattice interaction, and local Jahn-Teller distortions) are important in the orthorhombic manganites, while these effects are strongly reduced in the rhombohedral samples. The orthorhombic PI state was proposed to contain a superparamagnetic secondary phase in a certain temperature interval above $`T_c`$. These results are consistent with the x-ray diffraction measurements reported in our paper, providing further evidence of the correlation between the local structure and the average lattice symmetry, and thus revealing the physics underlying the drastic differences in the properties of manganites with different lattice symmetry.

We are grateful to K. I. Kugel, A. J. Millis, and J. Mira for important discussions. This work was supported by the NSF under grants No. DMR-0093143, and DMR-0103858, by the DOE under Contract No. DE-AC02-98CH10886, and by the NSF MRSEC program, Grant No. DMR-0080008.


format: text

text:
# PROPAGATION OF SECONDARY ANTIPROTONS AND COSMIC RAYS IN THE GALAXY

## INTRODUCTION

The result of our recent analysis (Moskalenko et al., 2001b, 2002), in agreement with calculations of other authors (Molnar and Simon, 2001), was that matching the secondary/primary nuclei ratio B/C using reacceleration models leads to values of the spatial diffusion coefficient apparently too large to produce the required $`\overline{p}`$ flux, when the propagated nucleon spectra are tuned to match the local $`p`$ and He flux measurements. Assuming the current heliospheric modulation models are approximately right, we have the following possibility (Moskalenko et al., 2003) to reconcile the B/C ratio with the measured flux of secondary $`\overline{p}`$’s. The spectra of primary nuclei as measured in the heliosphere may contain a fresh local “unprocessed” component at low energies thus decreasing the measured secondary to primary nuclei ratio. This component would have to be local in the sense of being specific to the solar neighbourhood, so that the well-known “Local Bubble” phenomenon is a natural candidate.

The low-density region around the Sun, filled with hot H i gas, is called the Local Bubble (hereafter LB) (e.g., Sfeir et al., 1999). The size of the region is about 200 pc, and it is likely that it was produced in a series of supernova (SN) explosions, with the last SN explosion occuring approximately 1–2 Myr ago, or 3 SN occuring within the last 5 Myr. Most probably its progenitor was an OB star association (e.g., Ma íz-Apellániz, 2001; Berghöfer and Breitschwerdt, 2002). Knie et al. (1999) report about an excess of 60Fe measured in a deep ocean core sample of ferromanganese crust suggesting a SN explosion about 5 Myr ago at 30 pc distance and the deposition of SN-produced iron on earth. Sonett et al. (1987) report an enhancement in the CR intensity dated about 40 kyr ago, which is interpreted as the passage across the solar system of the shock wave from a SN exploding about 0.1 Myr ago. It could also be that “fresh” LB contributions from continuous acceleration in the form of shock waves (Bykov and Fleishman, 1992), and/or energetic particles coming directly from SN remnants still influence the spectra and abundances of local CR.

## THE CALCULATION PROCEDURE

We use the propagation model GALPROP111GALPROP model including software and data sets is available at http://www.gamma.mpe-garching.mpg.de/$``$aws/aws.html as described elsewhere (Strong and Moskalenko, 1998; Moskalenko et al., 2002); for the present purpose the 2D cylindrically symmetrical option is sufficient. The spatial diffusion coefficient is taken as $`D_{xx}=\beta D_0(\rho /\rho _0)^\delta `$; the diffusion in momentum space and other details can be found in our earlier papers. The nucleon injection spectrum of the Galactic component was taken as a modified power-law in rigidity (Jones et al., 2001), $`dq(p)/d\rho \rho ^\gamma /\sqrt{1+(\rho /2)^2}`$, for the injected particle density. The heliospheric modulation is treated using the force-field approximation. The LB spectrum is taken in the form (Bykov and Fleishman, 1992, for continous acceleration by interstellar shocks): $`df/d\rho \rho ^\eta \mathrm{exp}(\rho /\rho _b)`$, where $`\rho `$ is the rigidity, and $`\rho _b`$ is the cut off rigidity parameter. In terms of kinetic energy per nucleon $`E`$ this can be re-written as $`df/dE=a(Z,A)\frac{A(E+m)}{Zp}\rho ^\eta e^{\rho /\rho _b},`$ where $`a(Z,A)`$ is the abundance of a nucleus $`(Z,A)`$, $`Z,A`$ are the nucleus charge and atomic number, correspondingly, $`m`$ is the atomic mass unit, $`p`$ is the momentum per nucleon, $`\rho _b=\frac{A}{Z}\sqrt{(E_b+m)^2m^2}`$. We show the results obtained with $`\eta =1`$, but the particular spectral shape of the LB component is not important as long as it decreases sharply towards high energies and is much softer than the Galactic CR spectrum.

The procedure to tune the CR elemental abundances, secondary/primary nuclei ratios, and $`\overline{p}`$ flux we adopted was as follows. The high energy part of B/C ratio *and* $`\overline{p}`$ flux measurements are used to restrict the parameters of the diffusion coefficient $`D_0`$ and $`\delta `$, while the low energy part of B/C ratio is used to fix a value for the reacceleration level and define the parameters of the LB component. The Galactic CR elemental source abundances are tuned using the CR abundances at high energies where the heliospheric modulation is weak. The LB isotopic abundances are tuned to match the low-energy data by ACE and Ulysses.

Figure 1 (left) illustrates the process of fixing the normalization of the diffusion coefficient using $`\overline{p}`$’s. The $`\overline{p}`$ flux shown is calculated with $`\delta =0.47`$ and different normalizations in the diffusion coefficient, $`D_{\mathrm{\hspace{0.17em}0}}=2.6,3.3,4.3\times 10^{28}`$ cm s-2 at $`\rho =3`$ GV (for $`\overline{p}`$’s this corresponds to kinetic energy $`2`$ GeV). The injection index $`\gamma `$ is taken equal to 2.28, and the Alfvén speed $`v_A=23`$ km s-1. The $`\overline{p}`$ flux at maximum, $`2`$ GeV, appears to be quite sensitive to the value of the diffusion coefficient and allows us to fix it at $`D_{\mathrm{\hspace{0.17em}0}}=(3.3\pm 0.8)\times 10^{28}`$ cm s-2. The exact value of $`\delta `$ is not critical since we compare with the $`\overline{p}`$ measurements at maximum. Inelastically scattered $`\overline{p}`$’s, the “tertiary” component, appears to be important at low energies only in the ISM. Figure 1 (right) shows corresponding calculations of the B/C ratio. A halo height of $`z_h=4`$ kpc is used (Strong and Moskalenko, 2001; Moskalenko et al., 2001a).

The local interstellar (LIS) $`p`$ and He spectra used in our $`\overline{p}`$ calculations are the best fits to the data as described in Moskalenko et al. (2002). Because of the measurements with large statistics (BESS – Sanuki et al. 2000, AMS – Alcaraz et al. 2000), and weak heliospheric modulation above 10 GeV, the error arising from uncertainties in the primary spectra is only $`5`$%. We note that as in the case of other nuclei, there should be an LB contribution to $`p`$ and He spectra. The Galactic injection spectra of $`p`$’s and He should thus be significantly flatter below several GeV to match the data points at low energies. This does not influence the $`\overline{p}`$ production because (i) the LB does not produce a significant amount of secondaries, and (ii) the $`\overline{p}`$ threshold production energy is high, $`10`$ GeV.

Further tuning can be done using the high energy part of the B/C ratio, which is not influenced by heliospheric modulation and supposedly contains only a Galactic component of CR. Index $`\delta 0.47`$ is chosen as giving the best match. A pure Kolmogorov spectrum, $`\delta =1/3`$, thus seems to be excluded by $`\overline{p}`$ spectrum data taken in combination with B/C ratio data. The low energy B/C data were used to fix the LB spectral parameter at $`E_b=500`$ MeV, but $`E_b=400600`$ MeV also provides good agreement with B/C data. In this way, we have been able to obtain a model simultaneously fitting $`p`$, He, $`\overline{p}`$, and B/C data.

## ABUNDANCES IN COSMIC RAYS AND COSMIC RAY SOURCES

The Galactic CR source elemental abundances are tuned (at a nominal reference energy of 100 GeV), by a least squares procedure, to the abundances measured by HEAO-3 (Engelmann et al., 1990) at 5.6, 7.5, 10.6, 16.2 GeV/n. Relative isotopic abundances at the source are taken equal to solar system abundances (Anders and Grevesse, 1989). The key point in the fitting procedure is to obtain the correct abundance of B.

The LB elemental abundances are tuned simultaneously with spectra using the low energy part of the B/C ratio and isotopic abundances at 200 MeV/n from ACE (Wiedenbeck et al., 2001) and at 185 MeV/n from Ulysses (DuVernois and Thayer, 1996). For many elements ACE and Ulysses abundances differ by 10%. For this reason, to the statistical errors we added 5% systematic error. The fitting procedure is influenced by the adopted value of the modulation potential. The following pairs of modulation potentials give almost the same final result: (ACE, Ulysses) = (400 MV, 700 MV), (450 MV, 600 MV), (500 MV, 500 MV). The data deviate from calculations in both directions, which mean that we are unlikely to introduce essential systematic error by assuming a wrong value of the modulation potential.

The diffusive reaccleration model with an LB component shows good overall agreement with data including secondary to primary ratios, spectra, and abundances. The derived Galactic CR source abundances and LB source abundances are plotted in Figure 2 relative to the solar system abundances (Grevesse and Sauval, 1998). The important result is that CR source and LB source abundances of all major elements (C, O, Ne, Mg, Si, S, Ar, Ca, Fe, Ni) are in good agreement with each other. Abundances of the Si group, Ca, Fe, Ni are near the solar system abundances. Abundances of other elements in Galactic CR and LB sources are mostly consistent with each other, and with solar system abundances, within a factor of 2. Relative to silicon, C, N, O, Ne, S are underabundant both in CR source and LB. This corresponds to the well-known first-ionization-potential or volatily correlation (see, e.g., Meyer et al., 1998).

Spectra of B, C, O, Fe are shown in Figure 3 for modulation levels 450 and 800 MV. In case of C, the normalization coefficient in the LB component is fixed as $`a(6,12)=6.35\times 10^4`$ cm-2 s-1 sr-1 for $`\eta =1`$. Since the elemental abundances are tuned at both high and low energies, the calculated sub-Fe/Fe ratio also agrees well with data.

Secondary nuclei K, Sc, Ti, V appear to be overabundant in the LB sources relative to the solar system (shown as upper limits in Figure 2), though the derived absolute LB abundances are not large. The derived LB abundance of Ti does not exceed that of Cr while the derived abundances of Sc and V are not larger than that of Mn. One possible reason for this excess is the uncertainty in the production cross sections, which is especially large for these nuclei. Sometimes there is no measurement at all, while phenomenological systematics are frequently wrong by a factor of two or more. Often, there is only one measurement at $`600`$ MeV/n, which has to be extrapolated in both directions. This allows only a nearly flat Webber-type or Silberberg-Tsao-type extrapolations while the real cross sections usually have large resonances below several hundred MeV, and decrease with energy above a few GeV (see, e.g., Moskalenko et al., 2001a). We note that Davis et al. (2000) used semiempirical cross sections based on Webber et al. (1990) and also predicted fluxes of sub-Fe elements which are too low.

An estimate of the overall error, which is reflected in the derived LB source abundances, can be obtained by assuming the absence of K, Sc, Ti, V in the LB source. Such a discrepancy is below 20%, and can be *removed* by allowing the production cross sections to increase at low energies by $`1520`$%, which seems plausible. Another possibility is errors in flux measurements of the rare CR species. Ulysses and ACE measurements are not always in agreement. Note that even for such an abundant nucleus as Fe, which is the main contributor to the sub-Fe group, the discrepancy exceeds 10%, while the disagreement in abundance of Sc is $`30`$%.

The derived source overabundance of sub-Fe elements in the LB could also in principle arise from composition differences between the ISM in the LB and solar or Galactic average ISM. This is suggested by the fact that the relative abundances of secondary elements in the LB sources are systematically larger than in the Galactic CR sources (Figure 2).

### Isotopic distributions in cosmic rays

Be and B isotopes are assumed all secondary, thus there is no possibility to tune them. Our calculation shows perfect agreement with the data on relative isotopic abundances of Be and B. This is in contrast with a “standard” reacceleration model, where we obtained a 15% discrepancy with relative abundances of Be and Be isotopes (Strong and Moskalenko, 2001). Abundances of stable isotopes of other elements are not very conclusive because they are present in the sources, but O and Si isotopic distributions still agree very well with data by ACE, Ulysses, and Voyager assuming only 16O and 28Si isotopes are present in the LB component. C and N isotopic distributions do not agree too well, but this may point to a problem with cross sections. The calculated ratio 13C/12C $`0.11`$ at 120 MeV/n ($`\mathrm{\Phi }=500`$ MV) in the model with LB contribution is still a factor $`1.5`$ too large compared to the measured value, $`0.0629\pm 0.0033`$ (Voyager 50–130 MeV/n, Webber et al., 1996) and $`0.078\pm 0.011`$ (Ulysses 100–200 MeV/n, DuVernois et al., 1996a), which may be connected in part with overproduction of 13C by 15N. If we replace the cumulative cross section 15N$`+p`$13C with cross section 14N$`+p`$13C, the calculated ratio 13C/12C will be lowered by 10% as estimated. Assuming the absence of the isotope 13C in the Galactic CR sources gives another 10% reduction. Altogether these corrections yield 13C/12C $`0.09`$, close to the data.

In case of radioactive isotopes Be, Al, Cl, usually used as “clocks” in CR, the agreement with the most accurate low energy data by ACE is very good and all the ratios are consistent with each other indicating that a halo height $`z_h=4`$ kpc is a good estimate. Higher energy data by ISOMAX (Be) are also consistent with calculations considering the large error bars. The 54Mn/Mn ratio indicates a somewhat smaller halo, but this may be related to uncertainty in its half-life and/or production cross section. The half-life of 54Mn against $`\beta ^{}`$ decay is the most uncertain among the four radioactive isotopes – it is only one which is not measured directly. Apart from the half-life, a possible source of errors can be production cross sections of Mn isotopes. Only the reaction natFe$`+p`$54Mn on a natural sample of Fe has been measured well enough (see compilation by Mashnik et al., 1998).

## CONCLUSION

In a previous paper we have shown that the $`\overline{p}`$ flux and B/C ratio appear to be inconsistent with measurements when computed in standard diffusion/reaccelation models. In this paper we have demonstrated that this discrepancy can be resolved if some part of the local CR consists of a “fresh” component accelerated in the LB. The independent evidence for SN activity in the solar vicinity supports this idea.

Combining the measurements of the $`\overline{p}`$ flux *and* B/C ratio to fix the diffusion coefficient, we have been able to construct a model consistent with measurements of important nuclei ratios in CR and derive elemental abundances in the LB. Calculated isotopic abundance distributions of Be and B are in perfect agreement with CR data. The abundances of three “radioactive clock” isotopes in CR 10Be, 26Al, 36Cl are all consistent and indicate a halo size $`z_h4`$ kpc based on the ACE data. 54Mn indicated a smaller halo, but this may be related to its half-life uncertainty and/or cross section errors. The derived fraction of LB component in CR is small compared to Galactic CR and has a steep spectrum with a cutoff above several hundred MeV/n. Other experimental data (except maybe the overabundance of Sc, Ti, V) do not contradict this hypothesis. The derived source overabundance of sub-Fe elements in the LB may be caused by trivial uncertainties in the production cross sections, or could in principle arise from composition and/or evolution differences between the ISM in the LB and solar or Galactic average ISM.

Recently there has appeared some indication (Huang et al., 2001) that the atmospheric contribution to the $`\overline{p}`$ flux measured in the upper atmosphere is underestimated by $``$30%. This means that the flux of $`\overline{p}`$’s in CR in reality may be *lower* at the top of the atmosphere by at least 25–30%. If the latter is true, the reacceleration model (even without LB) could still be the best one to describe propagation of nucleon species in the Galaxy.

## ACKNOWLEDGEMENTS

The authors are grateful to M. Wiedenbeck for providing the ACE isotopic abundances. I. V. M. is grateful to the Gamma Ray Group of the MPI für extraterrestrische Physik, where a part of this work has been done, for hospitality. I. V. M. and S. G. M. acknowledge partial support from a NASA ATP grant.

E-mail address of I.V. Moskalenko: imos@milkyway.gsfc.nasa.gov

Manuscript received ; revised ; accepted


text:
# Could the 55 Cancri Planetary System Really Be in the 3:1 Mean Motion Resonance?

## 1 Introduction

To date, 87 planetary systems, containing 101 giant extrasolar planets, have been discovered in Doppler surveys of solar-type stars (Butler et al. 2003), among which there are 10 multiple-planet systems, including 8 two-planet systems (HD 82943, GJ 876, HD 168443, HD 74156, 47 Uma, HD 37124, HD 38529 and HD 12661) and 2 three-planet systems (55 Cnc and $`\upsilon `$ And). As a second newly-discovered triple planet system, the most interesting characteristic of 55 Cnc is that the ratio of the orbital periods of the inner two companions (Marcy et al. 2002) is $`44.27/14.653`$, which suggests they are near 3:1 Mean Motion Resonance (MMR). However, the 3-Keplerian fits and other results given by Marcy et al. (2002) do not imply any evidence of the resonance for 55 Cnc at present. Here we examine revised orbital elements from Fischer et al.(2003) to search for evidence of a 3:1 resonance. Next, if the two planets are in the resonance, what kinds of the mechanisms make the system stable? In the present work, we aim to investigate the problem.

We completed pure N-body codes (Ji, Li & Liu 2002) to perform numerical simulations for 55 Cnc system by using RKF7(8)(Fehlberg 1968)and symplectic integrators(Feng 1986;Wisdom & Holman 1991). The time step is 0.36 d, about $``$2.5 percent of the orbital period of the innermost planet,which is sufficiently small for the integration. In addition, the numerical errors were effectively controlled with the accuracy of $`10^{14}`$ over the timescale from 1 Myr to several million years.

## 2 Dynamical Analysis

Two orbital solutions for the 55 Cnc have appeared in the literature: the first data set was reported by Marcy et al. (2002), and with more observations, Fischer et al. (2003) further updated the orbital parameters of 55 Cnc system. Their publications allow us to obtain the values for the integration. In this paper, we adopt the orbital parameters of the system provided by Fischer et al. (2003). In the simulation, the mass of the parent star is $`1.03M_{}`$, and those of three orbiting planets are adopted as $`0.83M_{jup}`$, $`0.20M_{jup}`$ and $`3.69M_{jup}`$ respectively (see Table 1), by assuming $`\mathrm{sin}i=1`$. In preparation of the orbital elements for integration, we assume that all of the planets are initially in the same plane, then fix their semi-major axes and eccentricities and the remainder angles (nodal longitude, argument of periastron, mean anomaly) are randomly generated for each orbit. In total, we carried out hundreds of the simulations for various initial conditions for the 55 Cnc system. More than 90 percents of the systems that are not linked to the 3:1 MMR become unstable in thousands of years. In Table 1, we list one set of the values out of 8 orbital solutions in association with the 3:1 MMR and apsidal resonance.

To examine our orbital solution based on Fischer et al. (2003), we used our orbits to compute the radial velocity. In Fig.1, the upper panel shows the comparison of the measured radial velocities and theoretically calculated radial velocities for 55 Cnc-the sign diamond represents the observations given by Marcy et al. (2002)(see their Table 2), and the sign plus denotes the results from our orbital solution; next we look at the lower panel which exhibits the residual velocity with rms of 11.50 m s<sup>-1</sup>.

### 2.1 3:1 Mean Motion Resonance

Many authors have studied the mean motion resonance in the multiple-planet systems: GJ 876 (Laughlin & Chambers 2001; Kinoshita & Nakai 2001; Lee & Peale 2002; Ji, Li & Liu 2002) and HD 82943 (Gozdziewski & Maciejewski 2001) are respectively in the 2:1 MMR, 47 Uma (Laughlin, Chambers & Fischer 2002) is close to a 7:3 commensurability. And we are now concentrating on the dynamics of the 3:1 resonance of 55 Cnc. The critical arguments for the 3:1 MMR are $`\theta _1=\lambda _13\lambda _2+2\stackrel{~}{\omega }_1`$, $`\theta _2=\lambda _13\lambda _2+2\stackrel{~}{\omega }_2`$ and $`\theta _3=\lambda _13\lambda _2+(\stackrel{~}{\omega }_1+\stackrel{~}{\omega }_2)`$, where $`\lambda _1`$, $`\lambda _2`$ are the mean longitudes of Companions B and C respectively, and $`\stackrel{~}{\omega }_1`$, $`\stackrel{~}{\omega }_2`$ respectively denote their periastron longitudes (hereafter subscript 1 denotes Companion B; 2, Companion C).

Marcy et al. (2002) pointed out that they did not find the evidence of the 3:1 MMR by showing the time behavior of $`\theta _1`$ with triple-Keplerian fit and self-consistent fits, for no libration of the resonant argument was detected. They also argue that further observations and dynamical fitting will clarify the relationship between the two inner planets of 55 Cnc. By contrast, we performed the integration on the basis of the fitting by Fischer et al. (2003). Figures 2a–2c exhibit the arguments of $`\theta _1`$, $`\theta _2`$ and $`\theta _3`$ plotted against time. From the figures, we see that $`\theta _1`$ librates about $`180^{}`$, $`\theta _3`$ walks around $`0^{}`$ with large amplitudes, and $`\theta _2`$ wanders near $`90^{}`$ with an amplitude of $`\pm 70^{}`$ for the timescale of 3 Myr. Because of the existence of the 3:1 resonance, as seen in Fig.3, the upper panel shows the semi-major axis of Companion B $`a_1`$ executes small vibrations about 0.115 AU, while the lower panel displays that of Companion C $`a_2`$ varies about 0.240 AU for millions of years. Obviously, these results confirm that the two companions of 55 Cnc are now really in the 3:1 MMR. Here we simply present the fact of the resonance in the sense of the numerical results, our future work will make an analytical theory to explain the orbital resonance for this planetary system.

However, one may ask what the situation was like before the two companions enter into the resonance. At first, we recall a recent work on the likely origin of GJ 876 by Lee & Peale (2002). They revealed the differential planet migration due to planet-nebular interaction could lead to capture into the 2:1 MMR. If we do not consider the outermost Companion D with a long period of 4780 days, the 55 Cnc system is quite similar to GJ 876, such as the semi-major axes and the masses of resonant pairs. In fact, the real system is not the case, because of the presence of Companion D, the process of the migration would possibly become more complicated. Future studies (Peale, in private communication) will help unearth the enigma of the origin of the resonance for 55 Cnc.

### 2.2 Secular Dynamics of 55 Cnc

For the system with small eccentricities and inclinations, the long-term variations of the eccentricities (or inclinations) of the planets can be described by the classical Laplace-Lagrange(L-L) secular theory provided that the mean motion resonance is absent. Now we draw attention to the secular behavior of 55 Cnc system. In Table 2 are given the data from the secular theory by considering all three planets. We see the contributions from Companion D is much less than one percent of the mutual effects of two inner planets. Hence, we can discuss the secular dynamics of Companions B and C by ignoring Companion D, without changing the quantitative results.

The linear L-L secular solutions for the eccentricities(Murray & Dermott 1999)can be written:

$$\begin{array}{c}e_j\mathrm{sin}\stackrel{~}{\omega }_j=e_{j1}\mathrm{sin}(g_1t+\beta _1)+e_{j2}\mathrm{sin}(g_2t+\beta _2),\hfill \\ e_j\mathrm{cos}\stackrel{~}{\omega }_j=e_{j1}\mathrm{cos}(g_1t+\beta _1)+e_{j2}\mathrm{cos}(g_2t+\beta _2),\hfill \end{array}$$
(1)
where $`e_{ji}`$, $`g_i`$ and $`\beta _i`$ ($`j,i=1,2`$, $`j=1`$ for Companion B, $`j=2`$ for Companion C) can be evaluated by the initial conditions in Table 1, and $`t`$ denotes time. For this case, we obtain two eigenfrequencies $`g_1=0.727^{}`$/yr and $`g_2=0.128^{}`$/yr. Together with $`\beta _1=4.92^{}`$, $`\beta _2=23.46^{}`$, and other determined constants in equation (1), we can derive the analytical expressions of the eccentricities $`e_j`$ for both of the inner planets. The numerical eccentricities of Companions B and C are respectively in the intervals of \[0.0, 0.22\] and \[0.20, 0.41\], while the theoretical values belong to those of \[0.0, 0.17\] and \[0.30, 0.41\],respectively. Because the inner two planets are now trapped in the 3:1 resonance, the secular theory is expected to be departure from the numerical results (we will construct a more detailed model under the consideration of the 3:1 MMR to study the system in our forthcoming paper), but it does correctly demonstrate the qualitative long-term variations. Therefore, we still perform integrations for the system that only consists of two inner planets to numerically examine the long-term orbital evolution without Companion D. We again find the outermost planet will not greatly influence on the secular motion of the inner companions because Companion D is fairly distant from Companions B and C(the ratio of $`a_2/a_30.04`$, $`a_2,`$$`a_3`$ are the semi-major axes of Companions C, D) and the mutual perturbations are very small. Furthermore, with additional, longer integrations of three planets, we find that the eccentricity of Companion D ranges from 0.26 to 0.29 and its secular motion is also steady. No sign shows this system is chaotic by numerical results as well as secular theory. Hence, we can safely conclude that the 55 Cnc system is dynamically stable for the lifetime of the star.

To explore the stability of 55 Cnc with respect to the planet masses, we fix $`\mathrm{sin}i`$ in increment of 0.1 from 0.1 to 0.9. In the complementary experiments, we integrate the system for 1 Myr with rescaled masses and the orbital parameters in Table 1. As a result, we find the system become unstable in short time when $`\mathrm{sin}i<0.5`$, indicating more massive planets are not fit for this system.

### 2.3 Apsidal Resonance

Apsidal resonance occurs when two planets share the same time-averaged rate of apsidal precession and the apsidal longitudes exhibit phase-locking. The relative apsidal longitudes is $`\mathrm{\Delta }\omega =\stackrel{~}{\omega }_1\stackrel{~}{\omega }_2`$. Malhotra (2002) described a dynamical mechanism for establishing apsidal alignment in a pair of planets that initially move on nearly circular orbits and showed that an impulse perturbation that imparts an eccentricity to one of the orbits excites the other planet’s eccentricity. Now it is widely believed that the apsidal resonance plays an important role in the stability of the discovered multiple-planet systems. In the case of $`\upsilon `$ And, the $`\mathrm{\Delta }\omega `$ of the outer two planets of the system librates about $`0^{}`$ (Chiang, Tabachnik, & Tremaine 2001), also in other systems-GJ 876 (Lee & Peale 2002) and 47 Uma (Laughlin, Chambers & Fischer 2002). The apsidal anti-alignment of HD 82943, librating around $`180^{}`$, was discovered by Ji & Kinoshita (in preparation) and it still occurs for HD 12661 (Lee & Peale, 2003). However, Beauge et al. (2002) argue that there may exist the planetary configurations of $`\mathrm{\Delta }\omega `$ librate neither $`0^{}`$ nor $`180^{}`$. In the case of 55 Cnc , we found that the $`\mathrm{\Delta }\omega `$ of Companions B and C librates about $`250^{}`$ with small amplitude for the timescale of 3 Myr (see Fig. 2d), and the phenomenon is never reported before and seems to be a paradigm of the statement by Beague et al. (2002). In comparison with the outer pairs of $`\upsilon `$ And, 55 Cnc differs in that the inner two planets in the apsidal resonance are simultaneously locked into the state of the 3:1 MMR. It is likely for the above two mechanisms to be responsible for the stability of the 55 Cnc system.

Next, what we care is whether the numerical libration of $`\mathrm{\Delta }\omega `$ of 55 Cnc could be explained in the analytical viewpoint. Thus we use the L-L secular theory to compute $`\stackrel{~}{\omega }_1`$ and $`\stackrel{~}{\omega }_2`$ in equation (1) respectively to show the predicted relative apsidal longitudes. In Fig. 4, we note the period of $`\mathrm{\Delta }\omega `$ from theory $`T_1=2\pi /(g_1g_2)600`$ yr and that of the numerical integration $`T_2`$ is almost 200 yr. The difference of two periods again results from the existence of 3:1 MMR. In addition, it is noteworthy that both of them have the libration width of the amplitude of $`\pm 60^{}`$. Moreover, by comparing two spans of libration, we can see the shift of the magnitude of $`70^{}`$ for the equilibriums, the minimums and maximums of the intervals. In summary, the analytical curve of the libration confirms the apsidal resonance found in numerical simulations.

## 3 Summary and Discussions

In this work, we obtain the orbital solutions of the planets of 55 Cnc in the 3:1 resonance. In the simulations, we find that not only three resonant arguments librate respectively, but $`\mathrm{\Delta }\omega `$ also librates about $`250^{}`$. These results imply the existence of the 3:1 MMR and the apsidal resonance for the studied system. We further emphasize that they can act as two important mechanisms of stabilizing 55 Cnc system. Moreover, it is the first time to observe the relative apsidal longitude librates neither $`0^{}`$ nor $`180^{}`$. We comment that the planetary systems in the mean motion resonance always accompany apsidal alignment, which can play a vital part in maintaining other non-resonant multiple-planet systems as well.

On the other hand, we argue that the libration of the resonant arguments may not always occur for even longer timescale, then the circulation temporarily takes place for short time and repeat the progress as time elapses, just as Marcy et al. (2002) pointed out the 3:1 resonance may appear in the past. However, the progenitor of the 3:1 MMR for 55 Cnc is not known today. Many works on the hydrodynamics simulations based on the interaction between the planets and disk (Bryden et al. 1999; Kley 2000; Papaloizou, Nelson& Masset 2001) are helpful to uncover the mystery.

We thank Stan Peale and Xinhao Liao for insightful discussions. We are grateful to Debra Fischer for helpful comments to improve the manuscript. This work is financially supported by the National Natural Science Foundations of China (Grant 10203005, 10173006, 10233020) and the Foundation of Minor Planets of Purple Mountain Observatory.


text:
# 𝐷-Meson Production from Recombination in Hadronic Collisions

## I Introduction

Since their discovery in the 1970’s charm quarks have proved valuable probes of Quantum Chromodynamics (QCD) via the production of open- and hidden-charm states in hadronic collisions at high energy . With the charm-quark mass, $`m_c1.5`$ GeV/$`c^2`$, being on the borderline of perturbative and nonperturbative momentum scales, mechanisms associated with both regimes are likely to play a role in their elementary production. In addition, nonperturbative effects are inevitably involved in subsequent hadronization, being related to the wave functions of the final-state mesons or baryons.

The usual way of describing charmed hadron production is based on the factorization theorem , i.e., a perturbative treatment for the creation of a $`c\overline{c}`$ pair through annihilation of partons within the incoming hadrons (folded over appropriate parton distribution functions), followed by individual hadronization of the $`c`$-quarks utilizing a (more or less) universal fragmentation function $`D(z)`$ (where $`z`$ is a normalized momentum fraction of the $`c`$-quark). Both stages in this description are not free of problems, or at least restricted to certain kinematic regimes. E.g., in leading-order (LO) perturbative QCD (pQCD), the elementary $`c\overline{c}`$ production cross section underestimates open charm production in hadronic collisions by large factors of $``$5 (the so-called $`K`$-factors). Next-to-leading order (NLO) estimates typically account for half of the deficit, with appreciable uncertainties inherent in underlying parameters such as the bare charm quark-mass, renormalization and factorization scales . This obviously leaves room for additional contributions of nonperturbative origin.

Concerning the hadronization of $`c`$ and $`\overline{c}`$ quarks, some empirical observations, especially at low transverse momentum ($`p_t`$), cannot be captured by standard fragmentation functions. E.g., flavor systematics at forward rapidities exhibit the so-called leading-particle effect ; that is, charmed hadrons containing light quarks corresponding to projectile valence quarks are much enhanced over their antiparticles. A natural explanation is provided by a recombination  of projectile valence with forward-produced (anti-) charm quarks , facilitated by small rapidity differences. More indirectly, it has been found  that $`x_F`$\- and $`p_t`$-distributions of $`D`$-mesons are well reproduced by delta-function fragmentation functions, i.e., the produced charm quarks do not seem to lose momentum as one would expect for a (realistic) fragmentation process. This is suggestive for the importance of recombination mechanisms even at central rapidities. As opposed to high-$`p_t`$ emission, charm quarks at low $`p_t`$ move at low velocity relative to their environment, with little energy available for string breaking, thus facilitating a ”statistical” attachment to neighboring light quarks.

Finally, as has been pointed out in Ref. , available experimental data from both $`p`$-$`N`$ and $`\pi `$-$`N`$ collisions point at a significant enhancement of $`D^+/D^{}`$ over $`D^0`$ total production (dominated by central rapidities and low $`p_t`$) over expectations from perturbation theory. This feature is also borne out when performing empirical fits using the PYTHIA event generator, which require separate $`K`$-factors to fit $`D^\pm `$ and $`D^0/\overline{D}^0`$, differing by a factor of $``$.

In the present article we address the afore-mentioned flavor dependencies by evaluating a suitably formulated recombination model for charged and neutral $`D`$-mesons extended to the central $`x_F`$ region. This allows to assess (”soft”) recombination processes involving sea-quarks from both projectile and target nucleons (or pions), and thus aims at a consistent treatment of forward and total flavor production asymmetries. In addition, one obtains estimates for $`D_s`$ meson production within the same framework.

The following presentation is organized as follows. In Sect. II we recall the main features of recombination mechanisms as have been applied in the literature before, and generalize it for our purpose of central production. LO pQCD $`c\overline{c}`$ cross sections are employed to obtain the primordial $`c`$-quark distribution, including a single (constant) $`K`$ factor adjusted to the measured total yields. The two key (nonperturbative) quantities to evaluate subsequent $`D`$-meson formation, namely a two-parton distribution function and a recombination function, are discussed in detail. The former provides the light quark for the ”coalescence” process, whereas the latter characterizes the pertinent overlap probability to form a $`D`$-meson. In Sect. III we compare our results with available data from both $`\pi `$-$`N`$ and $`p`$-$`N`$ reactions, and in Sect. IV we quote our predictions for upcoming measurements at RHIC in the $`p`$-$`p`$ mode at center-of-mass (CM) energy $`\sqrt{s}`$=200-500 GeV. We finish with conclusions and an outlook for future work in Sect. V.

## II Charm-quark recombination into $`D`$-mesons

### A General Outline of the Approach

The basic idea of the recombination approach is that quarks produced in hadronic collisions hadronize by ”coalescing” with (anti-) quarks pre-existing in the wave function of projectile (and/or target). This mechanism is quite different from the usual fragmentation, where produced quarks hadronize individually by string breaking, independent of their environment. Fragmentation functions therefore ought to be universal objects, which, in principle, can be extracted from $`e^+e^{}`$ collisions, where fragmentation is expected to be the sole source of hadron production. However, as already pointed out in the introduction, the application to hadronic collisions is limited. The first clear deviations have been identified in forward production of low-$`p_t`$ kaons and pions in $`p`$-$`p`$ collisions at ISR/Fermilab fixed target energies ($`E_{lab}100`$ GeV) , which showed strong enhancement of the $`\pi ^+/\pi ^{}`$ ratio at large $`x_F`$ (up to a factor $``$ 5), and even more pronounced for $`K^+/K^{}`$ . This lead Das and Hwa  to propose, within the parton model framework, the recombination model: forward produced quarks preferentially ”pick up” valence up-quarks from the projectile, thus favoring the ”leading” hadrons ($`\pi ^+=\overline{d}u`$ and $`K^+=\overline{s}u`$) over their ”non-leading” antiparticles. In addition, as shown in Ref. , recombination is also capable of nicely describing the observed increase of $`\pi ^{}/K^{}`$ and $`K^+/\pi ^+`$ for $`x_F1`$.

From a theoretical point of view, (anti-) charm quarks are a cleaner probe of recombination dynamics due to a negligible probability of producing them in secondary reactions. Indeed, the recombination approach has been successfully applied  (see also ) to forward production of charmed hadrons, which similarly exhibit leading particle effects. However, little attention has so far been paid to light-flavor asymmetries in the bulk production of charmed hadrons, which necessarily requires substantial contributions in the central region ($`x_F`$ around 0) at low $`p_t`$ where most of the yield is concentrated.

Here we generalize the recombination approach to incorporate $`D`$-meson formation at small $`|x_F|`$. As the low (Bjorken-) $`x`$ region is predominantly populated by sea-quarks, the main extension concerns the evaluation of their recombination with $`c`$-quarks. One of our ideas here is, that, besides the valence content, the proton sea possesses a well established flavor asymmetry which could reflect itself in the flavor composition of $`D`$-mesons.

To begin with, the evaluation of the production of $`c\overline{c}`$ quarks has to be specified. For forward $`D`$-meson production it has been suggested that an ”intrinsic” charm component  in the projectile wave functions could be responsible for significant recombination contributions . Such a ”hard” charm component seemed to be required to account for the leading $`D`$-meson distributions towards large $`x_F`$ at ISR energies. On the other hand, more recent data  indicate that the shape of inclusive $`D`$-meson $`x_F`$ distributions agrees well with pQCD predictions for bare (anti-) charm quarks (i.e., without momentum loss due fragmentation). We take this as a motivation for the following picture  that will be employed below: $`c\overline{c}`$ quarks are assumed to be exclusively created in (primordial) hard parton-parton collisions, evaluated in LO pQCD upscaled by a single empirical $`K`$ factor ($``$ 5-6 for typical choices of parton distribution functions (PDF’s) and charm-quark mass $`m_c=1.35`$ GeV) to match the experimental yields. The such generated charm-quark ($`x_F`$-) distributions are then subjected to recombination processes which necessarily involves nonperturbative information ***Note that in Ref. , a hard-scattered (projectile-) parton itself (after radiating off a gluon that fuses into $`c\overline{c}`$ with a target gluon) participates in subsequent recombination, whereas in the present analysis spectator anti-/quarks recombine. Formally, the former process is suppressed by one power in $`\alpha _s`$, which, for very forward production, is compensated by a large kinematic enhancement. Consequently, the results of Ref.  exhibit essentially no asymmetry for central $`x_F0.2`$.. For comparison with data, charm quarks that do not recombine are hadronized via isospin symmetric fragmentation with the usual polarization weights (i.e., 3:1 for $`D^{}/D`$), and, for simplicity, in $`\delta `$-function approximation for the $`x_F`$ dependence.

An important part of the inclusive charged and neutral pseudoscalar $`D`$-meson yields (as reported by most experiments) stems from feeddown contributions from $`D`$-meson resonances, most notably the vector-mesons $`D^{}(2010)`$. The two charge states of the latter exhibit rather asymmetric decay characteristics: the $`D^{}(2007)^0`$ decays to $``$ 100% in $`D^0(1865)`$ (with accompanying $`\pi ^0`$ or $`\gamma `$), whereas the $`D^{}(2010)^+`$ has a branching ratio of 68% into $`D^0(1865)\pi ^+`$ and only 32% into $`D^+(1869)`$ final states (charge conjugate modes for anti-$`D`$ mesons implied). As a consequence, the flavor composition of inclusive $`D`$-mesons depends on the relative abundance of directly produced $`D`$ and $`D^{}`$ mesons (higher resonances are neglected). For flavor-symmetric direct production, and defining $`ϵ=N_D^{}/(N_D^{}+N_D)`$, one finds for the inclusive charged-over-neutral $`D`$-meson ratio

$$R_{c/n}\frac{(D^++D^{})}{(D^0+\overline{D}^0)}=\frac{(1\frac{2}{3}ϵ)}{(1+\frac{2}{3}ϵ)},$$
(1)
whereas the inclusive particle-over-antiparticle ratios, $`R_{+/}`$ and $`R_{0/\overline{0}}`$, are independent of $`ϵ`$. For $`c`$-quark fragmentation, we employ the usual partitioning according to spin-isospin polarization counting implying a 3:1 ratio for $`D^{}:D`$, i.e., $`ϵ=3/4`$ and thus $`R_{c/n}=0.33`$. On the other hand, the recombination contribution is essentially a coalescence process and therefore expected to be sensitive to the mass of the product. Taking as a guideline effective thermal weights as extracted from light hadron production in $`pp`$ collisions , with $`T_H=`$ 170-175 MeV, results in a ratio $`D^{}/D3(m_D^{}/m_D)^{3/2}\mathrm{exp}[(m_D^{}m_D)/T_H]3/2`$, i.e., $`ϵ=0.6`$, which will be adopted below. In the absence of any flavor asymmetries this value for $`ϵ`$ implies $`R_{c/n}=0.43`$. The $`x_F`$-dependence of particle/antiparticle asymmetries is often displayed in terms of the so-called asymmetry function, defined via the production ratios as

$`A_{/+}(x_F)`$ $`=`$ $`{\displaystyle \frac{N_D^{}N_{D^+}}{N_D^{}+N_{D^+}}}={\displaystyle \frac{R_{/+}1}{R_{/+}+1}}`$ (2)
$`R_{/+}`$ $`=`$ $`{\displaystyle \frac{1+A_{/+}}{1A_{/+}}},`$ (3)
and likewise for $`\overline{D}^0`$ vs. $`D^0`$ (for baryon projectiles, anti-$`D`$ mesons are the ”leading” particles).

For the remainder of this section we focus on the recombination contribution to $`D`$-meson formation. In the above specified approach, the pertinent cross section in hadronic collisions at high energies takes the form

$$x^{}\frac{d\sigma ^{rc}}{dx_F}=\frac{dx_{\overline{q}}}{x_{\overline{q}}}\frac{dz}{z}\left(x_{\overline{q}}z^{}\frac{d^2\sigma }{dx_{\overline{q}}dz}\right)(x_q,z;x_F),$$
(4)
and a similar expression for $`\overline{D}`$’s upon replacing $`c\overline{c}`$ and $`\overline{q}q`$. Here, $`x^{}=E/(\sqrt{s}/2)`$ and $`x=p_l/(\sqrt{s}/2)`$ denote the $`D`$-meson energy and longitudinal momentum, respectively, normalized to the incoming proton (carrying CM energy and momentum $`\sqrt{s}/2`$), $`z^{}`$ and $`z`$ those of the charm-quark, and $`x_{\overline{q}}`$ the momentum fraction of the (projectile) light antiquark participating in the recombination (note that, for nucleons, the antiquark necessarily arises from the sea). The double differential cross section for the simultaneous production of $`c`$ and $`\overline{q}`$ quarks is given by

$`x_{\overline{q}}z^{}{\displaystyle \frac{d^2\sigma }{dx_{\overline{q}}dz}}={\displaystyle \underset{m_c^2}{\overset{W}{}}}𝑑m_{c,}^2{\displaystyle \underset{i,j}{}}{\displaystyle \underset{x_1^{min}}{\overset{x_1^{max}}{}}}𝑑x_1{\displaystyle \frac{x_{\overline{q}}x_1f_{\overline{q}i}^{(2)}(x_{\overline{q}},x_1)x_2f_j(x_2)}{(x_1z_+)}}{\displaystyle \frac{d\widehat{\sigma }_{ijc\overline{c}}}{d\widehat{t}}}`$ (5)
with kinematic boundaries $`x_1^{min}=z_+/(1z_{})`$, $`x_1^{max}=1x_{\overline{q}}`$ (defining $`z_\pm =\frac{1}{2}(z^{}\pm z)`$), and $`W=s(1zx_{\overline{q}})(1x_{\overline{q}})(1+z)/(2x_{\overline{q}})^2`$. $`d\widehat{\sigma }_{ij}/d\widehat{t}(\widehat{s},m_{c,})`$ denotes the standard LO pQCD parton fusion cross section  (including the $`K`$-factor). The summation over $`i,j`$ accounts for the possible parton combinations from target and projectile hadrons, with $`f_j`$ the distribution function for the target (nucleon). An analogous expression describes recombination processes with target (light-) quarks, upon replacing $`f_{\overline{q}i}^{(2)}(x_1)f_i(x_1)`$ and $`f_j(x_2)f_{\overline{q}j}^{(2)}(x_2)`$. The two crucial ingredients in evaluating the cross section, Eq. (4), are the recombination function $``$ and the two-parton distribution function $`f_{\overline{q}i}^{(2)}`$. In the following we discuss both entities in more detail.

### B Recombination Function

The recombination function $``$ essentially represents the (absolute value squared) wave function of the nascent charm-meson in terms of its quark constituents. For forward production one can rather straightforwardly employ the parton model in the collision (CM) according to

$$_D=\left[B(a,b)\right]^1\xi _{\overline{q}}^a\zeta ^b\delta (1\xi _{\overline{q}}\zeta )$$
(6)
with $`\xi _{\overline{q}}x_{\overline{q}}/x`$, $`\zeta z/x`$. The Beta-function prefactor, $`B(a,b)\mathrm{\Gamma }(a)\mathrm{\Gamma }(b)/\mathrm{\Gamma }(a+b)`$ ($`\mathrm{\Gamma }`$: Gamma-function), ensures the correct wave function normalization,

$$1=\underset{0}{\overset{1}{}}\frac{d\xi _{\overline{q}}}{\xi _{\overline{q}}}\underset{0}{\overset{1}{}}\frac{d\zeta }{\zeta }_𝒟.$$
(7)
The actual values of the exponents $`a`$ and $`b`$ have been motivated in two ways. On the one hand, in Ref. , the requirements that the average momentum fractions of the two constituents scales with their masses, and that $`G_D_D/(\xi _{\overline{q}}\zeta )`$ stays finite for $`\xi _{\overline{q}}0`$ (in the spirit of a constituent quark distribution function), are met by the choice $`a`$=1 and $`b`$=5. On the other hand, in Ref. , it has been argued on the basis of leading Regge trajectories that $`a=1\alpha _v`$ with $`\alpha _v`$=0.5 for $`v`$=$`u`$, $`d`$, and $`b=1\alpha _c`$ with $`\alpha _c`$$``$-2.2. Our results shown below turn out to have little sensitivity to the difference between these two cases.

A more severe complication arises when applying the parton model wave function for small $`D`$-meson momenta, i.e., for $`x_F`$ close to zero (central production); here, light-cone momentum fractions become ill-defined, involving even contributions from backward moving charm or light (anti-) quarks (w.r.t. the collision CM). Due to the, in principle, invariant nature of the parton-model we can resolve this problem by evaluating the recombination function in a boosted frame, which we choose as the projectile rest frame. In practice this procedure is, however, beset with a residual frame dependence owing to the finite transverse masses of charm- and, more sensitively, light-quarks. We will use the bare-charm quark mass together with typical light quark transverse masses of $`m_{q,t}0.30.5`$ GeV.

Another approach which avoids this complication is based on the formulation of the wave function in terms of rapidity variables, as first suggested in Ref.  via a (normalized) Gaussian distribution,

$$_D(y;y_c,y_q)=\frac{1}{\sqrt{2\pi }\sigma _y}\mathrm{exp}\left(\mathrm{\Delta }y^2/2\sigma _y^2\right).$$
(8)
The latter is solely characterized by its width $`\sigma _y`$ and the rapidity difference $`\mathrm{\Delta }yy_cy_{\overline{q}}`$ with the charm- and light-quark rapidities given by

$$y_i=\frac{1}{2}\mathrm{ln}\left(\frac{E_i+p_{i,z}}{E_ip_{i,z}}\right)$$
(9)
Again one is sensitive to the transverse mass, $`m_{q,t}`$ of the light quark within the incoming hadron, which here, however, carries a physical meaning: it is a nonperturbative quantity (indirectly) accessible by various experimental information (e.g., transverse momentum spectra of Drell-Yan pairs).

In Sect. III we will discuss both approaches of evaluating the recombination function, i.e., parton-model vs. rapidity-space wave functions. The pertinent results presented turn out to be rather robust. This, after all, can be understood from the notion that also the parton model wave functions imply the main contribution from a rather limited relative rapidity interval of $`\mathrm{\Delta }y`$=1-2 .

### C Two-Parton Distributions

#### 1 Kinematic Correlations

The two-parton distribution function (2-PDF) $`f_{\overline{q}i}^{(2)}`$, which importantly figures into the double-differential cross section, Eq. (5), characterizes the probability to simultaneously find parton $`i`$ at $`x_1`$ (participating in the hard fusion process into $`c\overline{c}`$) and the antiquark $`\overline{q}`$ with momentum fraction $`x_{\overline{q}}`$ (participating in the recombination process). In general it can be cast into the form

$$x_{\overline{q}}x_1f_{\overline{q}i}^{(2)}(x_{\overline{q}},x_1)=x_{\overline{q}}\widehat{f}_{\overline{q}}(x_{\overline{q}})x_1f_i(x_1)\rho (x_{\overline{q}},x_1),$$
(10)
where $`\widehat{f}_{\overline{q}}(x_{\overline{q}})`$ represents the constrained probability for parton $`\overline{q}`$ at momentum $`x_{\overline{q}}`$ under the condition that parton $`i`$ is at $`x_1`$, and $`\rho (x_{\overline{q}},x_1)`$ encodes all other correlations which are usually assumed to be governed by phase space.

As a simple ansatz the following factorized form has been proposed ,

$$f_{\overline{q}i}^{(2)}(x_{\overline{q}},x_1)=C_2f_{\overline{q}}(x_{\overline{q}})f_i(x_1)(1x_{\overline{q}}x_1)^p,$$
(11)
where $`C_2`$ is a normalization constant. The last factor in Eq. (11) implements a phase space dependence for large momentum fractions $`x_{\overline{q}}+x_11`$ (for a collinear 3-body phase where, in the CM system, target and forward $`[i\overline{q}]`$ state recede from each other with the projectile at rest, one has $`p=1`$ ). Some confidence in this ansatz was drawn from phenomenological successes in describing, e.g., (semi-) inclusive $`\pi `$\- and $`K`$-production in the fragmentation region . In subsequent work  the 2-PDF has been elaborated on a more microscopic basis within the so-called valon model of (low-$`p_t`$) hadronic structure, where light (anti-) quark distributions are generated from composite valence quarks at low momentum scales, see, e.g., Ref.  for a recent update. Within this approach, also $`c\overline{c}`$ production arises from Fock components of the valon distributions, similar to the ”intrinsic-charm” model . Their contribution can be significant at large $`x`$, but cannot account for inclusive charm yields. However, since our primary objective here is the evaluation of $`D`$-meson flavor dependencies in bulk production, we retain the $`c\overline{c}`$ cross section as a hard process in factorized form, cf. Eq. (5). Thus, results to be discussed below will refer to the simplified ansatz, Eq. (11), for the 2-PDF’s.

A more explicit treatment anti-charm quark recombination with leading (valence) light quarks within a factorized scheme has been undertaken in Refs. . At the price of more schematic underlying 1-PDF’s, a 2-PDF of the form

$$f_{vi}^{(2)}=Ax_v^{\alpha _v}x_1^\gamma (1x_1x_v)^m(1x_1)^k$$
(12)
has been proposed. Here, $`x_v`$ denotes the momentum fraction of the (recombining) valence quark in the nucleon, and the exponents $`\alpha _v`$, $`\gamma `$, $`k`$ and $`m`$ reflect the appropriate limiting behaviors for $`x_{1,v}0,1`$ (e.g., $`\alpha _v=0.5`$ for valence quarks, $`\gamma =1`$ for sea quarks and gluons). The ansatz, Eq. (12), is furthermore constrained by the normalization condition

$$\underset{0}{\overset{(1x_1)}{}}𝑑x_vf_{vi}^{(2)}(x_1,x_v)=f_i(x_1),$$
(13)
which can be used to fix the numerical constant $`A`$ as well as $`k`$ and $`m`$ according to the parton type $`i`$ (note that this can only be applied for the case of valence quarks $`v`$). In particular, it has been argued  that for fixed $`x_v`$, $`f_{vi}^{(2)}(\eta _1)`$ with $`\eta _1x_1/(1x_v)`$ can be interpreted as the single-parton distribution of the quark $`i`$ in the hadron remnant $`\stackrel{~}{h}_1`$ (with the valence quark $`v`$ removed), and that for $`x1`$ and $`zx`$ an approximate factorization into two single-parton distributions holds, i.e.,

$`f_{vi}^{(2)}`$ $``$ $`Cx_v^{\alpha _v}(1x_v)^{m\gamma }f_i(\eta _1)`$ (14)
$``$ $`f_v(x_v)(1x_v)^\gamma f_i(\eta _1).`$ (15)
In an alternative view, one might consider $`f_{vi}^{(2)}`$ at fixed $`x_1`$ (being determined on a shorter timescale – hard $`c\overline{c}`$ production – than $`x_v`$) as the single-parton distribution of quark $`v`$ in the hadron remnant without parton $`i`$. In this case, one straightforwardly obtains from Eq. (12) the factorized form

$$f_{vi}^{(2)}=\stackrel{~}{A}f_v(\eta _v)f_i(x_1)(1x_1)^1$$
(16)
with $`\eta _vx_v/(1x_1)`$ and $`\stackrel{~}{A}=1`$ upon imposing Eq. (13). This result suggests an immediate generalization for sea quarks according to

$$f_{\overline{q}i}^{(2)}=f_{\overline{q}}(\eta _{\overline{q}})f_i(x_1)(1x_1)^1$$
(17)
with $`\eta _{\overline{q}}x_{\overline{q}}/(1x_1)`$.

To illustrate the differences between the simple factorization ansatz with phase space correction, Eq. (11), and the newly suggested ansatz ”remnant” 2-PDF, Eq. (17), we show in Fig. 1 partice-antiparticle asymmetries in $`D`$-meson production for fixed target $`p`$-$`p`$ collisions. The results include the recombination contribution only, and are based on a rapidity-space recombination function, Eq.(8), with $`\sigma _y=0.5`$ (using GRV94-LO  for the underlying 1-PDF’s). In the central region, where both $`x_{\overline{q}}`$ and $`x_1`$ are small, mutual (anti-) correlations are suppressed, and the two ansätze indeed approximately coincide. However, the ”remnant”-2-PDF, Eq. (17), develops a somewhat steeper increase of the leading-particle asymmetry with increasing $`x_F`$. Note also, that recombination with valence quarks is the prevalent process even for central $`x_F`$ (indicated by nonzero values at $`x_F=0`$, as well as the larger asymmetry in the neutral ratio involving $`u`$-quarks). The $`\overline{d}/\overline{u}`$ asymmetry in the proton sea, which favors $`D^+`$ formation, does not significantly suppress the charged ratio. One should also point out that the absolute magnitude of the recombination cross section for the individual $`D`$-meson states is about 20% larger with the 2-PDF of Eq. (17) (which leads to larger total asymmetries once combined with isopsin-symmetric fragmentation contributions).

#### 2 Flavor Correlations

An obvious (global) flavor correlation that will be imposed throughout consists of ”vetoing” a light (anti-) quark flavor $`i`$ for the recombination process if this type has been used in the hard process of creating a $`c\overline{c}`$ pair, i.e., $`\overline{q}i`$. The significance of this anti-correlation mainly establishes itself for leading $`D`$-mesons ($`c\overline{q}`$) in forward production with meson beams , where the valence anti-quark provides a large part of the $`c\overline{c}`$ yield.

A more delicate issue is the one of local correlations, relating to anisotropies of the quark and gluon distributions within a hadron. To illustrate their possible effects we will investigate one rather extreme scenario, which can be motivated, e.g., by the importance of instantons in hadron structure. The latter provide a qualitative explanation for the observed $`\overline{d}/\overline{u}`$ excess in the proton  (cf. Ref.  for a recent overview), as follows. Due to the specific properties of the instanton-vertex, a valence $`u`$-quark coupling to it has to be accompanied by an in- and outgoing $`d`$-quark. The $`d`$-quark line can either be provided by another valence-quark (which generates a strong $`ud`$-diquark binding in the nucleon), or be closed off by a condensate insertion (which, after all, generates the constituent mass of the $`u`$-quark). The second case implies that $`u`$-valence quarks are essentially accompanied by a $`d\overline{d}`$ ”cloud”, and $`d`$-valence quarks by a $`u\overline{u}`$ one, thus generating a 2:1 asymmetry in $`\overline{d}/\overline{u}`$, which roughly corresponds to the experimental value in the $`x10^1`$ region. At the same time, one expects the gluon cloud to be concentrated around the valence quarks. Imposing such a flavor correlation on the 2-PDF (e.g., for gluon fusion, in 2/3 of the cases $`c`$ and $`\overline{c}`$ can only recombine with either $`u`$ valence or $`d`$ and $`\overline{d}`$ sea-quarks), will be referred to as ”local flavor correlations” below.

## III Comparison with Experiment

### A Total $`c\overline{c}`$ Yields and $`K`$-Factors

It is well-known that LO pQCD calculations for $`\overline{q}q,ggc\overline{c}`$, coupled with LO PDF’s for the colliding hadrons, underestimate total charm production by an appreciable magnitude characterized by a typical factor $`K`$ 5. The latter depends rather little on the specific parton distribution function, but bears some sensitivity to the value of the charm-quark mass. NLO results can in principle account for a substantial part of the discrepancy, but large uncertainties due to the choice of renormalization and factorization scales persist  (see also Ref.  for a recent update). We here adopt the LO results for $`d\widehat{\sigma }_{ijc\overline{c}}/d\widehat{t}`$ in Eq. (5), with the (bare) charm-quark mass fixed at $`m_c=1.35`$ GeV. For comparison we will allow for two different sets of PDF’s in our calculations, i.e., the GRV-94 LO  and MRST01 the MRST01 parametrization is performed to higher order (HO) in $`\alpha _s`$, which is, strictly speaking, not consistent with the LO evaluation of the $`c\overline{c}`$ cross section; on the other hand, nonperturbative effects introduced via the recombination mechanism evade any expansion in $`\alpha _s`$ for the $`D`$-meson production cross section of Eq. (5). sets for nucleons, as well as GRS-99 LO  for pions. Fig. 2 shows the fit to fixed target $`pN`$ (left panel) and $`\pi N`$ (right panel) reactions, which are reasonably well reproduced, albeit with significantly different $`K`$-factors (which might be related to the fact that for pion projectiles valence-antiquarks contribute appreciably through annihilation on target valence-quarks).

Based on the these total $`c\overline{c}`$ cross sections we will in the following address the magnitude and flavor composition of $`D`$-meson production from recombination based on Eqs. (4) and (5). As detailed in Sect. II, charm quarks that do not recombine will be hadronized with ”standard” fragmentation assuming no longitudinal momentum loss. We will mainly focus on $`x_F`$ distributions and excitation functions (i.e., inclusive yields for $`x_F>0`$ as a function of collision energy, $`\sqrt{s}`$) for the ratios $`(D^++D^{})/(D^0+\overline{D}^0)`$, $`D^{}/D^+`$ and $`\overline{D}^0/D^0`$.

### B Flavor Asymmetries I: Proton-Nucleon Collisions

#### 1 Inclusive $`x_F`$ Distributions

Let us first address the typical features of recombination cross sections as a function of $`x_F`$. Fig. 3 confronts E769 data  for inclusive $`D`$+$`D_s`$ mesons with our calculations for recombination with different wave functions $`_D`$ according to Eqs. (6), (8), and a common two-parton distribution function of factorized form, Eq.(11), with $`C_2^N=p_N=1`$ (at this point there is no sensitivity to more detailed properties of the 2-PDF). Also shown is (80% of) the pQCD yield for $`c`$ and $`\overline{c}`$ quarks (solid curve), assumed to convert into $`D`$-mesons without longitudinal momentum loss ($`\delta (xz)`$-fragmentation). This distribution (for the MRST01 set with $`K`$=5 as fixed above) reproduces the data reasonably well. In general, the fraction due to recombination is sizable in the central region (where it makes up 30-60% of the total cross section), and starts to dominate at forward $`x_F0.4`$, where it can even exceed the primordial $`c`$-quark distribution by ”pick-up” of a comoving light quark.

In more detail, we find that using the parton-model recombination function with $`a`$=1, $`b`$=5 and $`m_t=0.3`$ GeV gives a $``$50% contribution to the total $`D`$-meson yield. This is reduced to $``$30% upon using $`m_t=0.5`$ GeV (not shown), with most of the reduction occurring at central $`x_F0.2`$. On the contrary, with the rapidity space recombination function the sensitivity to the transverse mass is essentially absent, giving identical yields within 2% for the two values

of $`m_t`$. When increasing $`\sigma _y`$ from, e.g., 0.5 to 1, the $`x_F`$ distribution broadens slightly with the integrated yield again being stable (within 1%). Also, the use of rapidity wave functions increases the recombination fraction of $`c`$ and $`\overline{c}`$ quarks to almost 70%. Due to these rather robust features of $`_D(\mathrm{\Delta }y)`$, (and the absence of boost ambiguities for central $`x_F`$ as opposed to the parton model wave function , cf. Sect. II B) we will employ this form from now on (fixing $`\sigma _y=0.5`$) unless otherwise stated.

A similar comparison at higher energy is performed in Fig. 4 with nonstrange $`D`$-mesons from $`p`$-$`p`$ collisions measured by E743 . At $`\sqrt{s}`$=39 GeV, the fraction of the recombination cross section to the expected integrated yield (which corresponds to 2/3 of the total $`c\overline{c}`$ cross section) amounts to almost 80%, as compared to about 70 % for the corresponding calculation at $`\sqrt{s}`$=21.7 GeV, cf. Fig. 3. This is due to the fact that at higher energies smaller $`x`$ values of the parton distributions are probed where the occupancy is higher. Of course, for recombination fractions approaching 1 interference effects will play an increasingly important role, which have been neglected here. We also note that the recombination cross section into strange $`D_s`$-mesons, displayed by the dashed line in Fig. 4, constitutes about 25% of the nonstrange one, very similar to what has been inferred for total yields . This fraction is somewhat smaller when using the GRV94-LO distributions.

#### 2 Excitation Function of Flavor Ratios

From here on all our results for $`D`$-meson recombination cross sections are supplemented with (isospin-symmetric) $`\delta `$-function fragmentation of the remaining $`c`$ and $`\overline{c}`$ quarks (except at large $`x_F`$ where the recombination part may even exceed the $`c`$-quark distribution), with appropriate feeddown systematics from $`D^{}`$ resonances, as elaborated in Sect. II A. Throughout, the two parameters of the nucleon 2-PDF, Eq. (11), are fixed at $`C_2`$=0.8 and $`p`$=1, together with rapidity space recombination functions ($`\sigma _y=0.5`$). For the underlying 1-parton distribution functions we choose the GRV94-LO set (in order to be consistent with the GRS99  pion sets in subsequent sections). For nuclear targets, the appropriate neutron fraction has been employed.

Fig. 5 shows the pertinent $`\sqrt{s}`$-dependence of $`x_F`$-integrated nonstrange $`D`$-meson flavor ratios in $`pp`$ and $`pN`$ collisions, indicating reasonable agreement with available measurements for charged and neutral ratios (middle and lower panel, respectively), whereas the asymmetry in the charged over neutral ratio (upper panel) appears to be underpredicted. A marginal increase of the latter can be achieved by introducing local flavor correlations within a valence-quark type picture (as outlined in Sect. II C 2), but not enough to significantly improve the description of the data. For pure proton targets, an additional small increase for the charged over neutral ratio is found (dashed-dotted lines). In fact, for proton targets the increase in the individual recombination with both valence $`u`$-quarks to form $`\overline{D}^0`$ mesons and with sea $`\overline{d}`$-quarks to form $`D^+`$ mesons is very moderate, as indicated by the $`D^{}/D^+`$ and $`\overline{D}^0/D^0`$ ratios on proton targets (cf. dashed dotted lines in the middle and lower panels).

### C Flavor Asymmetries II: Pion-Nucleon Collisions

#### 1 $`x_F`$ Dependencies

Flavor asymmetries in forward $`D`$-meson production from $`\pi N`$ collisions are experimentally well-established . In Fig. 6 we compare results of our approach with the $`x_F`$-dependence of $`D^{}/D^+`$ and $`D^0/\overline{D}^0`$ ratios as measured by WA92  in 350 GeV $`\pi ^{}N`$ reactions.

The simple factorization for the 2-PDF, i.e., Eq. (11) with $`C_2=p=1`$ for pions, leads to a $`x_F`$-dependence of the charged ratio (upper dashed-dotted line) that is somewhat flatter than the data. The description is much improved (upper solid curve) when employing a weaker phase space suppression (reducing the exponent to, e.g., $`p_\pi `$=1/4) in connection with a smaller normalization constant (e.g., $`C_2^\pi =0.4`$). The latter suppresses the asymmetry in the central region, whereas the former entails a stronger increase towards forward $`x_F`$. On the other hand, the neutral $`D`$-meson ratio is insensitive to this modification. The main feature of these data is the near absence of a leading particle asymmetry, even at very forward $`x_F`$. The model calculations essentially reproduce this behavior, which arises from the fact that the $`\overline{u}`$-valence quark in the pion, which on average carries a large projectile momentum fraction, contributes significantly to the (hard) $`c\overline{c}`$ production process (via annihilation on $`u`$-valence quarks in the nucleon), and thus is not at disposal for subsequent recombination. For small $`x_F`$ the measured asymmetry seems to become even negative, which in the calculations is due to significant contributions from recombination processes of $`\overline{c}`$ quarks with $`u`$ quarks from target nucleons. This effect becomes more pronounced when introducing local flavor correlations within the nucleon (lower dashed line), which increases (decreases) recombination with $`u`$ ($`\overline{u}`$) into $`\overline{D}^0`$ ($`D^0`$) mesons; at the same time, the charged ratio becomes somewhat suppressed at negative $`x_F`$ (upper dashed line), due to reduced (enhanced) recombination with $`d`$ ($`\overline{d}`$) quarks into $`D^{}`$ ($`D^+`$).

At a higher projectile energy (500 GeV), the E791 collaboration measured the $`x_F`$-dependence of the inclusive $`D^{}/D^+`$ ratio. As before, the modified 2-PDF with $`C_2^\pi =0.4`$ and $`p_\pi =1/4`$ accounts much better for the data than with $`C_2^\pi =p_\pi =1`$, cf. Fig. 7. In general, over the range of fixed target energies discussed here, the projectile-energy variation in the $`x_F`$-dependencies of the flavor asymmetries is rather weak.

In Fig. 8 differential $`x_F`$-spectra for the absolutely normalized recombination + fragmentation model calculations (with $`K`$=3 in accordance with Fig. 2) for the individual $`D`$-meson states (underlying the solid lines in Figs. 6 and 7) are compared to data from WA92 . Reasonable agreement in all channels over the entire $`x_F`$-range is observed.

#### 2 $`\sqrt{s}`$ Dependencies

We finally turn to the excitation function of the flavor asymmetries in $`\pi N`$ collisions. Again, our combined fragmentation + recombination approach gives a rather satisfactory description of the observed ratios. In contrast to the $`pN`$ case (where, however, the data accuracy is rather limited), now also the charged over neutral ratio is in line with the data. Overall, the flavor asymmetries introduced by recombination processes improve the agreement as compared to contributions from (isospin-symmetric) fragmentation alone (represented by the dotted lines) in all channels.

It should also be noted that overall charm conservation has not been enforced explicitly. However, the slight excess of anti-charm quarks induced by the $`D^{}/D^+`$ asymmetry is conceivably compensated by a positive $`\mathrm{\Lambda }_c/\overline{\mathrm{\Lambda }}_c`$ asymmetry, which is also observed experimentally  (but which we did not address here). The same argument applies for proton projectiles, where both $`D^{}/D^+`$ and $`\overline{D}^0/D^0`$ are larger than one, but at the same time the $`\mathrm{\Lambda }_c`$ asymmetry is more pronounced than for pion beams .

We close this section with a note on $`D_s`$ mesons. For the recombination contribution alone, the ratio for (inclusive) strange over nonstrange $`D`$-mesons at $`\sqrt{s}=26`$ GeV turns out to be 12% for forward production (13% for all $`x_F`$), with almost no energy dependence. This ratio increases slightly to 14-15% when replacing the GRV-LO with the GRV-HO PDF’s. From a compilation in Ref. , the experimentally measured ratios can be found to be $`15.8\pm 2.9`$% for 230 GeV $`\pi ^{}`$ beams , $`16.7\pm 3.3`$% at 250 GeV , and $`11.6\pm 1.4`$% at 350 GeV, with a combined average of $`R_{D_s/D}=12.9\pm 1.2`$%. This value is in surprising agreement with the recombination approach. In fact, for the remaining charm quarks, fragmentation should apply. The prediction of the LUND event generator  (where hadronization is essentially performed by string fragmentation) turns out to be $``$8 %, somewhat below the data. However, the combined recombination + fragmentation result would still be close to the observed value.

## IV Predictions for $`p`$-$`p`$ at RHIC

An important part of the physics program at the Relativistic Heavy-Ion Collider (RHIC) at Brookhaven consists of (polarized) $`p`$-$`p`$ collisions up to a maximum CM of 500 GeV. Collisions at 200 GeV have been performed in the year-2 (2002) and year-3 running periods, although direct open charm measurements are not expected before detector upgrades have been performed. We here quote our predictions for $`D`$-meson flavor ratios by extrapolating our recombination + fragmentation approach according to Sect. III B to 200 GeV CM energy We recall that the recombination part of the charm cross section increases with collision energy; at 200 GeV, the calculated fraction of $`D`$+$`D_s`$-mesons amounts to 67% of the total $`c\overline{c}`$ yield (or 83% of the expected total $`D`$+$`D_s`$ yield; note that these numbers include the coefficient $`C_2=0.8`$ for the nucleon 2-PDF used in all results starting with Fig. 5); our estimates for the ratios at RHIC energy do not account for potentially sizable destructive interference effects.. We find:
$`R_{(D^++D^{})/(D^0+\overline{D}^0)}=0.40`$ ; $`R_{D^{}/D^+}=1.24`$ ; $`R_{\overline{D}^0/D^0}=1.35`$ ; $`R_{D_s/D}=0.23`$ .
These values are all rather close to the ones at fixed target energies.

## V Conclusions

Charm-meson production in hadronic collisions at fixed target energies has been investigated with emphasis on asymmetries in the light-flavor sector. The calculations have been performed within a recombination model for produced anti-/charm quarks with surrounding light quarks in projectile and target, supplemented with a standard fragmentation treatment for the remaining $`c`$ and $`\overline{c}`$ quarks. More precisely, the recombination component has been evaluated assuming a factorization of the perturbatively calculated (hard) $`c\overline{c}`$ production vertex and the subsequent (nonperturbative) ”coalescence” process. On the one hand, this requires the use of a 2-parton distribution function, for which we used a simple factorization ansatz (as first suggested 25 years ago). On the other hand, the recombination function, as compared to previous analyses, has been generalized to incorporate coalescence with sea-quarks within the colliding hadrons, necessary to address bulk production dominated by the central region ($`x_F0`$). Color correlations have been neglected throughout.

Employing recent parameterizations of parton distribution functions for nucleon target as well as pion and proton projectiles, we have found that (i) recombination processes contribute significantly to total $`D`$-meson yields and, (ii) the description of observed flavor asymmetries – $`(D^++D^{})/(D^0+\overline{D}^0)`$, $`D^{}/D^+`$ and $`\overline{D}^0/D^0`$ – in bulk production is improved throughout as compared to isospin-symmetric fragmentation alone. The nontrivial isospin asymmetries (i.e., beyond known feeddown corrections from $`D^{}`$ states) are driven by the valence- and sea-quark content of the colliding hadrons. No significant sensitivity was found with respect to local flavor correlations, as one might expect from valence-quark or quark-diquark substructures in the nucleon. Except for the charged-over-neutral asymmetry for proton projectiles (which is, however, beset with large uncertainties), the magnitude and (not very pronounced) collision energy dependence of available fixed target data is rather well reproduced.

As has been noted before, the same framework is also able to explain the differential $`x_F`$ dependencies of the flavor asymmetries as measured in $`\pi N`$ collisions. These become more pronounced at forward $`x_F`$, and are, within the model calculations, more sensitive to phase space correlations in the 2-parton distribution functions.

Future data on identified charm mesons at collider energies, such as expected from $`p`$-$`p`$ runs at RHIC, will be most valuable to further test the relevance of recombination mechanisms in the hadronization process, and thus contribute to a better understanding of the latter.

ACKNOWLEDGMENTS

We gratefully acknowledge useful discussions with G. Sterman. This work was supported in part by the U.S. Department of Energy under Grant No. DE-FG02-88ER40388.