A&A 394, 1039-1056 (2002)
DOI: 10.1051/0004-6361:20021176
D. Maurin ^{1} - R. Taillet ^{1,2} - F. Donato ^{3}
1 - Laboratoire de Physique Théorique LAPTH, 74941 Annecy-le-Vieux, France
2 - Université de Savoie, 73011 Chambéry, France
3 - Università degli Studi di Torino and INFN,
Torino, Italy
Received 18 June 2002 / Accepted 29 July 2002
Abstract
In a previous study (Maurin et al. 2001), we explored
the set of parameters describing diffusive propagation of cosmic rays
(galactic convection, reacceleration, halo thickness, spectral index and
normalization of the diffusion coefficient), and we identified those
giving a good fit to the measured B/C ratio.
This study is now extended to take into
account a sixth free parameter, namely the spectral index of sources.
We use an updated version of our code where the reacceleration
term comes from standard minimal reacceleration models.
The goal of this paper is to present a general view of the evolution
of the goodness of fit to B/C data with the propagation parameters.
In particular, we find that, unlike the well accepted
picture, and in accordance with
our previous study, a Kolmogorov-like power spectrum for diffusion
is strongly disfavored. Rather, the analysis points towards
along with source spectra
index 2.0.
Two distinct energy dependences are used for the source spectra:
the usual power-law in rigidity and a law modified
at low energy, the second choice being only slightly preferred.
We also show that the results are not much affected by a different
choice for the diffusion scheme.
Finally, we compare our findings to recent works, using other propagation models.
This study will be further refined in a
companion paper, focusing on the fluxes of cosmic ray nuclei.
Key words: ISM: cosmic rays
(1) |
(2) |
Most analyses of cosmic ray nuclei data assume given power-laws for the diffusion and acceleration energy dependence, so that the results partially reflect certain theoretical a priori. In this work, we try to avoid this bias by determining the quantities and directly from the data, in particular B/C, for reasons exposed below.
The paper is organized as follows. We first recall the main features of our diffusion model. As a few modifications have been made since previous works, Sect. 3 is devoted to their description and justification. Then, the analysis method is described in Sect. 4 and the results are shown and discussed in Sect. 5; a comparison is eventually made with other similar works in Sect. 6.
This paper and its companion (Donato et al. in preparation) use the same description of cosmic ray propagation as our previous analyses (Maurin et al. 2001; Donato et al. 2001; Donato et al. 2002; Barrau et al. 2002; Maurin et al. 2002). Particles are accelerated in a thin galactic disk, from which they diffuse in a larger volume. When they cross the disk, they may interact with interstellar matter, which leads to nuclear reactions (spallations) - changing their elemental and isotopic composition - and to energy losses. Interaction with Alfvén waves in the disk also leads to diffusive reacceleration. The reader is referred to Maurin et al. (2001) - hereafter Paper I - for all details, i.e. geometry and solutions of our two zone/three-dimensional diffusive model, nuclear parameters (nuclear grid and cross sections), energy losses terms (adiabatic, ionization and Coulomb losses), solar modulation scheme (force-field), as well as general description of the procedure involved in our fits to data (selection of a set of parameters, test comparison to data). In particular, though some inputs are modified (see next section), the final equation describing cosmic ray equilibrium is formally equivalent to that of Paper I (see Eq. (A13)): it is a second order differential equation in energy solved with the Crank-Nicholson approach (see Donato et al. 2001, Appendix B - hereafter Paper II). Finally, a schematic view of our diffusion model is presented in Barrau et al. (2002) and Fig. 1 (see next section) summarizes the algorithm of our propagation code.
Some aspects of this model are formally unrealistic. First, the distribution of interstellar matter has a very simple structure: it does not take into account a possible z distribution inside the disk (thin disk approximation is used instead), nor radial and angular dependence in the galactic plane. The orthoradial dependence would even be more important from an accurate description of the magnetic fields and the ensuing diffusion, as flux tubes are likely to be present along the spiral arms. However, this is not crucial as we are interested in effective quantities (diffusion coefficient and interstellar density) but not in giving them a "microscopic'' explanation. This is why we chose to use a universal form of the diffusion coefficient, with the same value in the whole Galaxy. Finally, it is known that a fully realistic model has to take into account interactions between cosmic ray pressure, gas and magnetic pressure, i.e. magnetohydrodynamics.
The semi-analytical diffusion approach should be thought of as an intermediate step between leaky box approaches and magnetohydrodynamics simulations and is actually justified by these two very approaches: the first showed that the local abundances of charged nuclei can be roughly described by two phenomenological coefficients - the escape length and the interstellar gas density in the box. The second hints at the fact that the propagation models such as the one used here are well suited for the description of cosmic ray physics.
However, it is difficult to conclude whether these parameters are valid for other kinds of cosmic rays (, , nuclei induced -ray production) and whether they are either meaningful but valid only locally on a few kpc scale (i.e. not in the whole Galaxy - see as an illustration Breitschwerdt et al. 2002), or meaningless but phenomenologically valid as an average description of more subtle phenomena (see as an example the discussion of the Alfvénic speed in Sect. 6.3.4).
(3) | |||
This equation can be rewritten using the cosmic ray differential
density
.
As the momentum distribution function is normalized
to the total cosmic ray number density (
),
we have
to finally obtain
From a theoretical point of view, the most natural choice for the energy dependence of the source term seems to be a power-law in rigidity (or momentum) for . This translates into in our set of equations (see Eqs. (4) and (5) above). Several different forms were used in the past because of the lack of strong evidence from observed spectra (see for example Engelmann et al. 1985; Engelmann et al. 1990). In particular, our previous analysis allowed only a rigidity dependence (for the special case ). These two forms differ only at low energy and we chose to keep them both to estimate their effect on our results. As we show below, it is quite small.
Finally, different diffusion schemes lead to different forms for the energy dependence of the diffusion coefficient and the reacceleration term. Several aspects of the diffusion process are treated in Schlickeiser (2002), and we considered three alternative possibilities: (i) Slab Alfven wave turbulence, with and , (ii) Isotropic fast magnetosonic wave turbulence, with and , and (iii) mixture of the two last cases, and .
All results will be presented with the case (i), except in the specific discussion in Sect. 5.5.
(8) |
Figure 1: Diagrammatic representation of the various steps of the propagation code. | |
Open with DEXTER |
This procedure is very time consuming. Even when the location of minima in the six-dimensional parameter space are known, more than configurations are needed to have a good sampling of the regions of interest, for a given form of the source term energy dependence.
(12) |
We also set the halo thickness L to 6 kpc, leaving us with
four free parameters (,
K_{0},
and ).
All curves depicted in Fig. 2 correspond to one-dimensional
cuts through the absolute
minimum (for a given ,
the
three different cuts justify the fact that we are located in a minimum).
In the upper panel of Fig. 2, we plot the values of the as a function of K_{0}/L, for different values of
(and the corresponding
).
Figure 2: Evolution of the value for various combination of parameters. All curves are for type (a) spectra ( ) with and the halo size is fixed to L=6 kpc. Each curve shows one-dimensional cuts in the K_{0}/L (upper panel), (left lower panel) or (right lower panel) direction of the 3-dimensional hyper-surface. In the upper panel is varied from 1.0 to 0.3 and, in the lower panels the same symbols to indicate are conserved. Each curve gives the absolute minimum for the parameter on the abscissa axis, L being fixed to 6 kpc (similar curves with slightly different minima are obtained for other L values). | |
Open with DEXTER |
In the lower panels we present two cuts in the two other directions, namely in the and directions. The first one tells us that except for the special case for which the curve skips to null , B/C is fitted with between 10 and 20 km s^{-1}. The best are for convective velocity around 16-18 km s^{-1}. For km s^{-1} (and ) the goodness of the fit quickly decreases. We can see that when is around 0.4-0.3, the B/C ratio becomes very sensitive to the values. It appears that when is decreased, a good fit is maintained provided that is also lowered. This is possible down to for which the best value for is zero. For lower , the previous trade-off cannot be achieved (as must be positive for the galactic wind to be directed outwards) and no good fit is possible.
The right panel shows the curves as functions of the Alfvén velocity. The minimization procedure always yields a far different from zero. Good fits are obtained for values of 40-50 km s^{-1}.
In each of the explored directions, the curves are very narrow: the diffusion model leads to meaningful and interpretable values for all the physical, free parameters. Similar results, with slightly different values for the minima, are obtained for the other values of L in the range kpc.
Figure 3: Left panel: evolution of the best value with K_{0}/L for various (1.0 to 0.3, from left to right) at a fixed . Each curve correponds to a given halo size L from 14 kpc to 2 kpc. Right panel: the same best values are presented versus and . In both panels, empty circles correspond to type (a) spectra and stars to type (b) spectra. | |
Open with DEXTER |
In Fig. 3 we present the results for the same analysis for different values of the halo thickness L and considering also the form (b) for the source spectra, i.e. . The total spectral index is still set to 2.8. The left panel reports the as functions of K_{0}/L, for different values of and L, and for both types of source spectra. We see that the choice (b) globally improves the fit, and the favoured range for is now (whereas for choice (a)). At fixed and L, the absolute minima for both choices correspond to very similar values of K_{0}/L. We can also notice that type (a) spectra are, for the higher , more sensitive to variations of L.
In the right panels we show a cut in the - plane. For both type (a) and (b) spectra, yields a null value for the convective wind. Type (a) spectra give a little bit higher . At fixed , the variation of L has almost no effect on , while it is strongly correlated with the increase of .
In this section we discuss the results obtained when the index is varied between 1.3 and 2.5, being set to a given value which has been extensively used in the literature.
Figure 4: Same as in Fig. 2 (type (a) spectra, L=6 kpc), but for a fixed . | |
Open with DEXTER |
Figure 4 corresponds to the previous Fig. 2. In the left panel we observe that a large variation of the index has a slight effect on the normalization of the diffusion coefficient K_{0}, which stays around an average value kpc Myr^{-1} for L=6 kpc. Evolution of the absolute minimum is also far less sensitive to than (see previous section). However, for the fit to the data is poor and a global power at is excluded.
The lower panels represent a cut in the and directions. We can observe that the minimization procedure always drives the minima towards convective velocities between 12 and 16 km s^{-1}, the least being obtained for the smallest . This range is again very narrow. Similarly, reacceleration is needed to fit data and the minima of the are obtained for between 55 and 75 km s^{-1}. Towards this lower limit, is high and the model cannot confidently reproduce observations.
When is fixed, we can conclude that a variation in the power of the type (a) source spectrum does not strongly act on the evolution of and also . This can be also easily understood: forgetting for a while energy gains and losses, we see from diffusion equation solutions (the same behavior occurs in leaky box models) that the source term can be factorized so that secondary to primary ratios finally do not depend on Q(E), i.e. are independent of . Once again, the absolute minimum is identified by a steep in these three directions.
Figure 5: Same as in Fig. 3 but for a fixed . | |
Open with DEXTER |
In Fig. 5 we present the results for , and for both type (a) and (b) source spectra, to focus on the evolution of L and . The left panel tells us that the evolution of the halo thickness from 2 to 14 kpc, at fixed (in other words, at fixed ) does not change the goodness of the fit. Only a slight modification in K_{0}/L is required in order to recover the same B/C flux ratio. Type (b) source spectra reproduce quite well the data for all the explored parameter space. On the contrary, the better theoretically motivated type (a) spectra cannot reproduce observations for if . Since at high energies the two source spectra are equivalent, we must conclude that it is the low energy part of B/C which is responsible for such a discrimination.
The right panels show the absolute minima in the - plane. Both spectra require non-null reacceleration and convection. Even more so, the selected values reside in the narrow interval for , i.e. -15 km s^{-1} and between 40 and 90 km s^{-1} for the Alfvén velocity.
Figure 6: Same as in Fig. 2 (type (a) spectra, L=6 kpc), but for a fixed . | |
Open with DEXTER |
Figure 6 describes the results of the analysis done assuming type (a) source spectra, with fixed index and L=6 kpc. A consensus seems to emerge in favor of values (see Drury et al. 2001 and references therein), close to the index given by primeval acceleration models, but any other value would be fine for the purpose of this section. In the upper panel has been varied between 1.0 and 0.3, and the figure shows the evolution of the with respect to K_{0}/L. As in Fig. 2 and, at variance with Fig. 4, the minima correspond to K_{0}/L spanning over almost two orders of magnitude. It is the modification of the power-law in the diffusion coefficient - and not in the source spectrum - that significantly acts on K_{0}. Once again, the Kolmogorov spectrum is disfavoured: in this case it is obvious that the calculated flux ratio would be too hard. The best fits are obtained for -0.9.
The lower panels show the cuts in the and directions. The left one tells us that for smaller , the preferred convective velocities are smaller (and the best is larger), down to for which a no-convection model is prefered, with a bad . The best fits are obtained for around 15-18 km s^{-1}. In the right panel we can notice, again, that only models with reacceleration have been chosen by the minimization procedure. Lower point to higher K_{0}/L and values and lower . The same trend is recovered in the other cases treated above. Reacceleration and convection act, in a certain sense, in competition, even if data always give preference to a combined effect rather than their absence.
This trend (the smaller ,
the larger K_{0}, or equivalently K_{0}/L
as L is constant in the above figures) was already mentioned in
Sect. 5.1. Actually, as we will see in Sect. 6, the correlation
between K_{0}/L and
is more properly explained by virtue of
Eq. (17) so that the evolution of
is fixed by the evolution
of the two other free parameters, i.e. K_{0}/L and .
As regards ,
it only appears in Eq. (9). A rough estimation
can be inferred using power-laws
and
in Eqs. (4) and (9):
Figure 7: Same as in Fig. 3 but for a fixed . | |
Open with DEXTER |
Figure 8: Best values for various L (2, 6 and 10 kpc) in the plane . Left histograms are type (a) spectra and right histograms type (b). Notice that for right histograms, only the upper figure displays the values and . They have been omitted in the two remaining figures to gain contrast (for any L, these configurations have ). Assuming L=6 kpc, type (a) source spectra give a best value for and whereas type (b) gives for and . These were obtained with 26 data points. | |
Open with DEXTER |
Figure 9: From top to bottom: for each best in the plane (L=6 kpc), the corresponding values of , and are plotted for both source spectrum types. | |
Open with DEXTER |
Figure 10: Best values, in the plane , for the three different forms of the diffusion coefficient and reacceleration terms (i) Slab Alfven wave turbulence, with and , (ii) Isotropic fast magnetosonic wave turbulence, with and , and (iii) mixture of the two last cases, and . | |
Open with DEXTER |
We know present the result of the full analysis, in which all the parameters are varied. Figure 8 shows the evolution of the in the and plane for different values of L. We can see that, at fixed type (a) or (b) spectra, a change in the halo height L has almost no effect on the best surface. Generally, high values for are preferred and, a Kolmogorov regime for the spatial diffusion coefficient is strongly disfavoured over all the parameter space. More precisely, type (b) spectra point towards a band defined by in the - plane, whereas the type (a) spectra gives the additional constraint (see Fig. 8).
In Fig. 9 we show the preferred values of the three remaining diffusion parameters K_{0}, and , for each best in the - plane, when L has been fixed to 6 kpc. The two upper panels show that the evolution of does not affect K_{0}. On the other hand, as already noticed, we clearly see the (anti)correlation between the two parameters K_{0} and entering the diffusion coefficient formula, giving the same normalization at high energy ( ). Almost the same numbers are obtained for type (a) and (b) spectra. K_{0} spans between 0.003 and 0.1 kpc^{2} Myr^{-1}. We will discuss in the following sections how these results can be compared to the literature.
The middle panels show the values for the convective velocity. Only very few configurations include , always when , for both types of source spectra. The value of increases with . For type (a) spectra, increasing and at the same time makes change its trend. As remarked previously, the effect of Galactic wind is more subtle since it acts at intermediate energies and is correlated with all the other diffusion parameters through the numerous terms of the diffusion equation.
The lowest two panels show the influence of . We recover a correlation similar to the one discussed for K_{0} (see Eq. (13)). The Alfvén velocity doubles from to 0.3, whereas it is almost unchanged by a variation in the parameter (or equivalently ).
All the three analysed parameters (i.e. K_{0}, and ) behave very similarly with respect to a change in the source spectrum from type (a) to type (b). It can be explained as the influence on the primary and secondary fluxes can be factored out (see Sect. 5.2) if energy changes are discarded (their effect is actually small on the derived parameters). Existing data on B/C do not allow us to discriminate clearly between these two shapes for the acceleration spectrum. This goal could be reached by means of better data not only for B/C but also for primary nuclei (Donato et al., in preparation).
As discussed in Sect. 3.1, we tested three different diffusion schemes, with three different forms for the diffusion coefficient. Most results are basically insensitive to the choice of this form. In particular, the figures corresponding to Fig. 9 are almost identical to the case presented above, so that they will not be reproduced here. Figure 10 displays the as a function of and . The values of are slightly different in the three cases, but the general trend is the same, and all the previous conclusions still apply.
In an ideal situation in which we had very good and consistent data on B/C and sub-Fe/Fe ratios, the best attitude would be to make a statistical analysis of the combined set of data. Unfortunately, this is not currently the case. We consider two ways to extract information from the Sub-Fe/Fe data. First, as a check, we compare the sub-Fe/Fe ratio predicted by our model - using the parameters derived from our above B/C analysis - with data from the same experiment. Second, we search directly the minimum of the sub-Fe/Fe ratio, with no prior coming from B/C. As previously emphasized (see Sect. 4.2), this procedure is more hazardous since the statistical significance of the sub-Fe/Fe data is far from clear.
For each set of diffusion parameters giving a good fit to the observed B/C ratio, the sub-Fe/Fe ratio can be computed and compared to the values measured by HEAO-3. This is not as straightforward as in the B/C case because although Sc, Ti and V - that enter in the sub-Fe group (as combined in data here) - are pure secondaries, some of the species intermediate between sub-Fe and Fe, contributing to the sub-Fe flux, are mixed species (i.e. Cr, Mn). As a consequence, all the primary contributions were adjusted so as to reproduce the sub-Fe/Fe ratio at 3.35 GeV/nuc. The sub-Fe/Fe spectra are not steep enough at high energy, so that normalization at 10.6 GeV (i.e. as for B/C) would have led to less good fits. We emphasize that to perform this normalization of secondary-to-primary is equivalent to making an assumption about the elemental composition of the sources, which is usually deduced from secondary-to-primary ratios. A different choice would slightly shift the normalization of sub-Fe/Fe ratio without affecting much our conclusions.
Figure 11 displays the values obtained when the diffusion parameters giving a good fit to B/C are used to compute the sub-Fe/Fe ratio, for each value of and (for type (a) spectra and L=6 kpc, although the results for type (b) and/or different L are quite similar). This surface is very similar to the surface obtained with B/C, pointing towards high values of (compare to Fig. 8).
We now consider a full sub-Fe/Fe analysis (i.e. the parameters minimizing are looked for) but we emphasize that the results given here are from our point of view far less robust than those obtained from B/C. As a consequence, conclusions of this section have to be taken only as possible trends. Several points can be underlined from Fig. 12: (i) as for the B/C case, the best is obtained for type (b) spectra. (ii) the general behavior of K_{0}, and to a less extent is mostly the same as for B/C. (iii) the type (b) spectra yield propagation parameters which are closer to B/C's, as we can see from values; (iv) finally, consistency with B/C analysis would be better obtained for pointing towards 0.6-0.7.
Typical spectra (modulated at MV) are shown in Fig. 13, for different values of the parameters and , along with the data points from HEAO-3 (Engelmann et al. 1990) and balloon flights (Dwyer & Meyer 1987). Three low-energy data points, from HET on Ulysses (Duvernois & Thayer 1996), HKH on ISEE-3 (Leske 1993) and Voyager (Webber et al. 2002) are also shown; they all have about the same modulation parameter, i.e. MV. The ACE points ( MV) are also displayed (Davis et al. 2002).
Figure 11: Values of obtained by applying the diffusion parameters - type (a) source spectra and L=6 kpc - giving the best fit to B/C, for each and , to sub-Fe/Fe. The are computed with HEAO-3 data points. | |
Open with DEXTER |
All the models displayed give similar spectra, which would be difficult to sort by eye. This may explain why some of these models (e.g. those with ) are retained in other studies. The main features are (i) the influence of on the high energy behaviour - a good discrimination between these models would be provided by precise measurements around 100 GeV/nuc - and (ii) the type (a) source spectra are steeper than type (b) at low energy.
To compare the reacceleration terms employed, we retain only
the spallation term and the highest order derivative in energy in
the diffusion equation, giving
Figure 12: From top to bottom: best and for each best in the plane (L=6 kpc), the corresponding values of , and are plotted for both source spectrum types. | |
Open with DEXTER |
Figure 13: The B/C and sub-Fe/Fe spectra (modulated at MV) for several sets of parameters (giving the best fit to B/C for these values) are displayed, along with experimental data from HEAO-3 (Engelmann et al. 1990), ballon flights (Dwyer & Meyer 1987), HET on Ulysses (Duvernois & Thayer 1996), HKH on ISEE-3 (Leske 1993) and Voyager (Webber et al. 2002). Note that ACE data (Davis et al. 2002) correspond to a modulation parameter MV. | |
Open with DEXTER |
Once this
rescaling - that differs from one
paper to another - is taken into account, a comparison
is possible between models if a minimal resemblance exists between the other input
parameters,
i.e. same ,
(plus same form of the source spectrum) and
halo size L; Table 1 shows the value adopted for these parameters
in two recent studies.
Maurin et al./This work | Seo & Ptuskin/Jones et al. | Moskalenko et al. | |
(2001)/(2002) | (1994)/(2001) | (2002) | |
Thin disk h (pc) | pc | h=200 pc | Gas distribution |
Halo size L (kpc) | -- | L=3 kpc | L=4 kpc |
Surface mass density | 10^{-3} g cm^{-2} | / g cm^{-2} | g cm^{-2} |
/ |
Jones et al. use the same set of equation and parameters
as Seo and Ptuskin for the stochastic reacceleration model. Consequently,
it seems that
defined in Jones et al. is the half-height of the
reacceleration zone, contrarily to what is depicted in their Fig. 4.
The surface mass density is defined as
where
is
the matter density in the thin disk.
The results are expected to be slightly different from our previous
study as the components have been modified.
First,
has a different interpretation in the two studies
(see Table 1, first column).
As underlined above - remembering that in Paper I the diffusion coefficients scaled
as
-, the Alfvén speed value
from Paper I (
)
has to be rescaled into
(i.e. as the standard convention
used in this work and others) through the relation
(15) |
However, we find that the conclusions raised in Paper I, in particular the behaviors reflected in Figs. 7 and 8 of Paper I, are basically unchanged (it is not straightforward to compare with present figures, but the careful reader can check this result using the above scaling relation and the corresponding parameter combinations). To be more precise, it appears that K_{0}/L does not significantly change (for example, for and L=2 kpc, we still have kpc Myr^{-1}, see Fig. 3 left panel - this paper - and Fig. 7 of Paper I). As regards the galactic convective wind, is shifted towards higher values, whereas the range remains roughly unchanged.
This can be easily understood: the additional term - comparable to a first order gain in energy, see Eq. (7) - has to be balanced to keep the fit good. This balance is ensured by enhanced adiabatic losses, i.e. bigger . Other parameters are only very slightly affected by this new balance.
Let us make a few comments at the qualitative level. First, starting with Mos02's models, we can note that convection has always been disfavored by these authors. For example, in their first paper of a series (Strong & Moskalenko 1998), a gradient of convection greater than 7 km s^{-1} kpc^{-1} was excluded. We notice that this result was not very convincing since it is clear from an examination of their figures that none of the models they proposed gave good fits to B/C data. Thanks to many updates in their code, their fits were greatly improved (Strong & Moskalenko 2001; Moskalenko et al. 2002) but if it is now qualitatively good, it is hard to say how good it is since no quantitative criterion is furnished. Anyway, our Fig. 2 allows us to understand why convection is disfavored in such models. Actually, if - as in the Kolmogorov diffusion slope hypothesis -, we see that for such a configuration, the best fits are obtained for km s^{-1}.
Similar comments apply to Jon01's models. Given a Kolmogorov spectral index for the diffusion coefficient, their combined fit to B/C plus sub-Fe/Fe data is not entirely satisfactory. It improves for higher values of and in the convective model (they do not include reacceleration in this model), their best fit being obtained for . As the authors emphasized, the search in parameter space was not automated and they cannot guarantee that their best fit is the absolute best fit. Actually, the sub-Fe/Fe contribution to the value has to be taken with care. First, the error bars are not estimated well enough to give a statistical meaning for values (see Sect. 4.2) and a different weight should be considered for B/C and sub-Fe/Fe. Second, if the best parameters extracted from B/C data reproduce formally the same surface when applied to the evaluation of sub-Fe/Fe (see Fig. 11), the direct search for the parameters minimizing for the same sub-Fe/Fe data gives constraints that are much weaker (see Fig. 12). Thus, any conclusion including this ratio is from our point of view far less robust.
Tables 2 and 3 give the results of Mos02 and
Jon01 - without any rescaling of any parameters -
compared to what is obtained here; only a few models are displayed.
L | h | K_{0} | Ref. | |||||
(kpc) | (kpc) | (g cm^{-2}) | (kpc) | (kpc^{2} Myr^{-1}) | (km s^{-1}) | (km s^{-1}) | ||
4. | n(r) | 1.6 | 4.0 | 0.201 | 0. | 30. | Good | Mos02 |
3. | 0.2 | 2.4 | 1.0 | 0.196 | 0. | 40. | 1.8 | Jon01 |
3. | 0.1 | 1.0 | 0.1 | 0.0535 | 0. | 105.8 | 4.2 | (Figs. 8 and 9, this paper) |
3. | 0.2 | 2.4 | 1.0 | 0.127 | 0. | 47.3 | 4.4 | This work |
For this model, the exact values are
,
.
This model give best fit to flux using a slightly modified
form for the source;
.
Corresponds to the best fit for the presented L,
and
value.
Table 2 shows K_{0}, and for and . Taking the first three lines at face value, our values of K_{0} and are very different from the others and our model seems to have a problem. However, the matter disk properties (height and surface density) are different in these models. To be able to compare, we set these quantities to the values given in Jon01 and the resulting parameters are shown in the last line of Table 2.
Actually, we know that in diffusion models,
the behavior is driven by the location of the closer edge, leading to a preferred
escape on this side. With L=3 kpc,
our three-dimensional model should behave as the two-dimensional
model with infinite extension in the r direction of Jon01.
This hypothesis can be validated if one takes their Eq. (3.6).
For the pure diffusion model (reacceleration and convection are discarded), one has a
simple relation between ,
L and K_{0} through an equivalent leaky box
grammage
A similar expression may be obtained in the presence of galactic wind:
in the wind model (their Eq. (4.6)), one has
L | h | K_{0} | Ref. | |||||
(kpc) | (kpc) | (g cm^{-2}) | (kpc) | (kpc^{2} Myr^{-1}) | (km s^{-1}) | (km s^{-1}) | ||
3. | 0.2 | 2.4 | 1.0 | 0.024 | 29. | 0. | 1.5 | Jon01 |
3. | 0.1 | 1.0 | 0.1 | 0.0056 | 15.5 | 35.3 | 3.0 | This work |
3. | 0.2 | 2.4 | 0.1 | 0.0134 | 36.5 | 26.5 | 3.1 | - |
Corresponds to the best fit for the presented L,
and
value.
Even with this rescaling, the diffusion coefficients obtained by the different authors quoted above are still not fully compatible. Another possible effect, namely the spatial distribution of cosmic ray sources, is now investigated. We note that in our model, the radial distribution of sources q(r) follows the distribution of supernovæ and pulsar remnants. The choice of this distribution has an effect on B/C spectra and on the parameters giving the best fits. If we use a constant source distribution q(r)=cte with set to 0 to follow Jon01, we find that K_{0} is enhanced by about 10%. We checked that it is also the case for results presented in Table 2. Hence, it appears that results for of Jon01, though slightly different, are not in conflict with ours. As regards and Mos2, using the scaling relation (16) along with a 10% decrease of K_{0} for Jon01, we obtain respectively K_{0}=0.226 (Mos02), 0.176 (Jon01), 0.127 (this paper) kpc^{2} Myr^{-1}. Thus there is some difference between Mos02 and Jon01, which is not obvious when the values taken naively from Table 2 are compared. These discrepencies could have several origins: treatment of cross-sections (we checked that total and spallative cross sections - taking into account ghost nuclei, see Paper I - are compatible with recent data, e.g. Korejwo et al. 2002), average surface density in Mos02 that is probably not exactly 1.6, choice of data and fits for Jon01 that differ from ours (some of the point they used are significantly lower than HEAO-3's). Finally, the fact that we scan the whole parameter space can make a difference from manual search. To conclude, results are qualitatively similar, but a few quantitative differences remain. The intrinsic complications and subtleties of the various propagation codes make it difficult to go further in the analysis of these differences.
The normalization K_{0} gives a measure of the efficiency of the diffusion process at a given energy. Its value can be predicted if (i) a good modelling of charged particles in a stochastic magnetic field and (ii) a good description of the actual spatial structure of this magnetic field, were available. It is not the case and the precise value of K_{0} is of little interest. Moreover, the presence of effects other than pure diffusion can be mimicked, at least to some extent, by a change in K_{0}. Equation (17) gives a whole class of parameters giving the same results and can be used to extract an effective value of K_{0} taking into account the effect of the size of the halo L and wind . This also explains the great range of values that can be found in the literature.
This relation shows that there is also an indeterminacy of the absolute density of the model, because as long as is constant, the grammage is also constant. Fortunately, a realistic distribution of gas can be deduced by more direct observational methods, so that a definite value of can be used.
We note that in our model, Galactic wind is perpendicular to the disk plane and is constant with z. Actually, the exact form of galactic winds is not known. From a self-consistent analytical description including magnetohydrodynamic calculations of the galactic wind flow, cosmic-ray pressure and the thermal gas in a rotating galaxy, Ptuskin et al. (1997) (see also references therein) find a wind increasing linearly with z up to kpc, with a z=0 value of about 22.5 km s^{-1}. Following a completely different approach, Soutoul & Ptuskin (2001) extract the velocity form able to reproduce data from a one-dimensional diffusion/convection model. They obtain a decrease from 35 km s^{-1} to 12 km s^{-1} for z ranging from 40 pc to 1 kpc followed by an increase to 20 km s^{-1} at about 3 kpc. For reference, our values for the best fits correspond to about 15 km s^{-1}. The difficulty to compare constant wind values to other z-dependences is related to the fact that cosmic rays do not spend the same amount of time at all z, so that there cannot be a simple correspondence (see also next section) from one model to another. As a result, all the above-mentioned models are formally different, with different inputs (spectral index, diffusion slope). Nevertheless, their values are roughly compatible, Ptuskin et al's model providing the grounds for a physical motivation for this wind. However, an even more complicated form of the Galactic wind could be relevant for a global description of the Galaxy (see Breitschwerdt et al. 2002).
(19) |
Second, the total rate of reacceleration (at least in a first approximation, see discussion below) is given by a convolution of the time spend in the reacceleration zone and the corresponding true Alfvén speed in this zone. There is a direct analogy with the case of spallations and the determination of the true density in the disk, as discussed above. The problem is still somewhat different, as there are no direct observational clues about the size of the reacceleration zone, or said differently, about . This leads to a degeneracy in that holds as far as , due to the structure of the equations in the thin disk approximation. For example, a model such as Strong et al's that uses cannot be simply scaled to ours. A cosmic ray undergoing reacceleration at a certain height z has a finite probability of escaping before it reaches Earth, this probability being greater for greater z. As a result, the total reacceleration undergone by a cosmic ray is actually not a simple convolution of the reacceleration zone times the Alfvén speed in this zone, but rather should be an average along z taking into account the above-mentioned probability (in principle, this remark also holds for the gaseous disk, though the latter is known to be very thin, a few hundreds of pc).
To conclude, there are basically three steps associated with three levels of approximations to go from the deduced from cosmic ray analysis to the physical quantity. First, if is approximatively constant with z, how large is the reacceleration zone height? The second level is related to the possibility that strongly depends on z in a large reacceleration zone. If it is too large, the link with the phenomenologically equivalent quantity in a thin zone is related to the vertical occupation of cosmic rays. However, this latter possibility seems to be unfavoured by MHD simulations (see Ptuskin et al. 1997). Finally, with the above parallel between interpretation of and , we see how misleading it is to obtain precise physical quantities from our simple model, since there is no one-to-one correspondence between reality and simplified models. This discussion shows that even if the actual derived Alfvén speeds are consistent with what is expected from "direct'' observation (10-30 km s^{-1}), the cosmic ray studies would certainly be not very helpful in providing physical quantities better than a factor of two. If we reverse the reasoning and retain our best models with L=6 kpc, we could conclude that must be 4 in order to give realistic values for (with evident a priori about ).
Starting with the leaky box; it has been shown more than thirty years ago (Jones 1970) that the concept of "leakage-lifetime'' was appropriate for the charged nuclei considered here (see also Jones et al. 1989), even if it broke down for (all orders in the development in "leakage'' eigenmodes contribute because of synchrotron or inverse Compton losses) and for radioactive nuclei (Prishchep & Ptuskin 1975). The leaky box, due to its simplicity, is very well suited for the extraction of source abundances (elemental as well as isotopic). It can also be used for secondary antiproton production, since the same processes as for secondary stable nuclei are at work. However, as emphasized in Paper II, leaky box models are not able to predict any primary contribution in the antiproton signal, since it requires the knowledge of the spatial distribution of primary progenitors. Considering a possible extension of leaky box models for stable charged nuclei to high energy (PeV), it has been demonstrated in Maurin et al. (2002) that they are to a good approximation sufficient to describe the evolution of cosmic rays. Last, it is well known that leaky box parameters are just phenomenological with only a distant connection to physical quantities.
This was further realized by Jones (1978, 1979) who first remarked that the phenomenological behavior of the escape length at low energy could be due to the presence of a Galactic wind. Jones et al. (2001) investigated further this idea and generated several equivalent phenomenological escape lengths from several possible physical configurations of a one-dimensional diffusion model. The relation between one-dimensional models and leaky box models is thus firmly established and very well understood. Moreover, this relation elucidates some of the physical content of leaky box models. Now if one wishes to overcome the inherent limitations of these models and say, to compute some primary antiproton component, one has to go through a three-dimensional model. It is likely that these models can also be related to the Jones et al. models (see Taillet & Maurin, in preparation). Several arguments used in the previous sections illustrate this view, but this occurs at least if the halo size is small compared to radial extension of the Galaxy.
In the semi-analytical two-zone model used here, it is possible to evaluate the primary antiproton component (see Barrau et al. 2001) and to take into account radioactive species, even in the presence of a local very underdense bubble (see Donato et al. 2002 for details). Our model fails to consider species such as and , since the latter suffer from large energetic losses in the halo so that no simple semi-analytical approach can be used. The parameters extracted from these models are much easier to interpret in terms of physical quantities.
Most of the limitations mentionned above are overcome by Strong et al's models. In this fully numerical model, all cosmic ray species can be computed self-consistently with the same propagation parameters. The main difference with our model is that a more realistic matter distribution is used instead of a thin homogeneous disk. They also consider that reacceleration occurs in the whole diffusion halo, which in our opinion is an approximation no more justified than the fact to confine it in a thin disk (see discussion in Sect. 6.3.4). Considering the gas distribution, both models are equally predictive for the charged nuclei (including antiprotons, see Fig. 9 of Paper II). On the one hand, our approach is better suited to scan the whole parameter space as we did in Paper I and in this paper. On the other hand, the Strong et al. models can check the consistency of , and with observations, and can include whatever deviation from ideal cases for K_{0}, , and more generally for any ingredient that enters in the description of propagation models.
To conclude about models and their use, Jones et al.'s approach is probably the best and simpliest way to understand how physical parameters affect the propagated flux. Our model is very well suited for a consistent evaluation of all charged nuclei and extraction of propagation parameters; furthermore it is an intermediate step where general behaviors can still be analytically explored (Taillet & Maurin, in preparation). In the Strong et al. model, all fine effects can be studied and modelled, with the counterpart that the numerical approach makes the physical intuition of the results less straightforward. In its present form, Strong et al's model can be viewed as a fully numerical version of ours, so that their behaviors are very close. This discussion could leave the reader with a feeling that apart from these different modellings left to personal taste, galactic propagation phenomena are well understood. It is surely not the case! Even if all these models are equivalent to describe the local observations of charged cosmic rays, they lead to very different conclusions and interpretations when the spatial variation of the cosmic ray density is considered. As an illustration of the poor current understanding of this global aspect, we mention the ever-lasting problem of the gamma ray excess about 1 GeV towards the Galactic center or the too flat radial -ray distribution observed in the disk (see Breitschwerdt et al. 2002).
For the rest, the conclusions of this paper can be summarized as follows: (i) we performed for the first time a full analysis of diffusion/convection/reacceleration models in the whole 6-dimensional parameter space (, , K_{0}, L, , ), and the values -0.9 and are preferred; (ii) this preference holds whatever the specific form of the spectrum at low energy; the numerical values of the other parameters are also only slightly modified by this low energy dependence even though deviation from a power-law at low energy is preferred. The study of fluxes should give a more definite answer; (iii) K_{0} scales logarithmically with and models with small halos tend to one-dimensional models with a simple relation between , K_{0}, L and (see also Taillet & Maurin, in preparation); (iv) several existing models are compared and the qualitative and quantitative differences between them are studied and partially explained.
Acknowledgements
D.M. would like to thank Aimé Soutoul for valuable remarks and Michel Cassé and Elisabeth Vangioni-Flam for many interesting discussions. We thank Prof. Schlickeiser for providing us with the diffusion coefficients in the different schemes in a useful form.