An Alternative Explanation of the Varying Boron-to-carbon Ratio in Galactic Cosmic Rays

Acceleration by strong shocks (compression ratio of 4 or more) is theoretically expected to produce ${E}^{-2}$ spectra during most of the life of a supernova remnant (Ellison & Eichler 1985), yet the Galactic cosmic ray spectral index is ${E}^{-2.7}$. Some of the discrepancy may be attributed to the fact that the last stages of shock acceleration are with shocks that have a compression ratio that is less than 4. Mostly, however, the discrepancy is attributed to an energy-dependent escape ratio, as is evidenced by the secondary-to-primary ratio, which decreases with energy. But the recent announcement (Aguilar et al. 2016) of an antiproton spectrum that is identical to the primary spectrum challenges this, because antiprotons¹ are also secondaries. Moreover, the positrons, which may also be secondaries, also have an identical spectrum to that of the primaries. Cowsik & Madziwa-Nussinov (2016) have recently suggested that the escape rate from the Galaxy is energy-independent, and that the primary source spectral index is really 2.7 and not 2. However, it would then be unclear why the boron-to-carbon (B/C) ratio is energy dependent. To address this question, Cowsik and coworkers (Cowsik & Wilson 1973) had already suggested that most of the boron is made in "nested" leaky boxes that encompass the production sites of cosmic rays and had also suggested that the escape rate from the nested leaky box, though not from the Galaxy, has the required energy dependence, but it is unclear why only one of the two escape rates would be energy dependent.

In this paper, I suggest an alternative reason for the decrease of the B/C ratio with energy, which does not demand an energy dependent escape rate. While the idea is somewhat speculative, I believe it should not be overlooked, even if it is only partly true, because it may play a role in determining the values of the measured quantities.

The Galactic disk is very thin. Most of the baryonic matter is concentrated within layers that are of the order of 100 pc from the equatorial plane, and this is only about 10⁻² of the Galactic radius R. The limits of CR anisotropy at $E\gtrsim 1$ TeV suggest that their streaming velocity is of the order of ${10}^{-4}c$ or less, meaning that the CRs do not stray far from their sources before escaping the disk (Parker 1969). See also D'Angelo et al. (2016) for an analysis of diffusion in self-generated turbulence that is consistent with this conclusion. Parizot & Drury (2000) have argued that cosmic ray secondary production near the production site of the primary elements can explain the evolution of Li, Be, and B in the Galaxy. The disk is also very inhomogeneous and different regions could have different densities, ionization fractions, different rates of star formation, different levels of Alfvén wave turbulence, etc. The correlation lengths of these varying quantities could be smaller or larger than the disk thickness. Accordingly, the cosmic ray escape rate from the disk, and the grammage of interstellar matter that they traverse before doing so, could differ among various parts of the disk.

The maximum energy to which a supernova shock can accelerate cosmic rays can also depend strongly on the Galactic environment where the shock exists. Ion-neutral damping, for example, can severely limit E_max (Bell 1978). The fact that cosmic rays appear to be made in regions where there are grains suggests that warm regions of the interstellar medium are well represented in overall cosmic ray production. Gamma rays from supernovae remnants display a wide variety of spectral indexes and cutoffs, illustrating the diversity of cosmic ray spectra that could be expected.

In this paper it is suggested that E_max at the site of cosmic ray production is anticorrelated with the grammage traversed by a cosmic ray before it escapes the Galaxy (or the subregion where most of the grammage is traversed). This would create an energy dependence in the B/C ratio even if the escape rate from the Galaxy is energy independent. The suggestion is quite general, but a specific version of it will be suggested following a more general description. We do not claim that the escape rate of cosmic rays from the Galaxy is entirely energy independent, but merely that an energy-dependent escape rate is not the only reason for an energy dependence in the B/C ratio.

Let us express the primary CR production spectrum at a given site of production as ${N}_{p}(E,{\boldsymbol{x}})$ below E_max, where ${N}_{p}(E)\,={E}^{-2-p},\,p\geqslant 0$, and where ${N}_{p}(E,{\boldsymbol{x}})$ vanishes at $E\geqslant {E}_{\max }$. The maximum energy E_max specific to that source is determined by any one of several physical considerations that will be discussed later. The overall spectrum ${N}_{p,T}(E,{\boldsymbol{x}})$ at any point ${\boldsymbol{x}}$ in the Galaxy is then ${N}_{p}S({E}_{\max }\geqslant E)$, where $S({E}_{\max }\geqslant E)$ is the total source contribution in sources that allow acceleration up to or beyond energy E. This fact can be written as

$$\begin{eqnarray}&&{N}_{p,T}(E,{\boldsymbol{x}})={E}^{-2-p}{\int }_{E}^{\infty }[{dS}({E}_{\max },{\boldsymbol{x}}){/{dE}}_{\max }]{{dE}}_{\max },\end{eqnarray} \tag{ 1 }$$

where ${dS}({E}_{\max },{\boldsymbol{x}}){/{dE}}_{\max }$ is the relative strength of sources with a maximum energy between E_max and ${E}_{\max }+{{dE}}_{\max }$ contributing to point ${\boldsymbol{x}}$. Note that the "relative strength" $S({E}_{\max },{\boldsymbol{x}})$ depends both on the number and strength of regions with maximum energy E_max as well as their distance from point ${\boldsymbol{x}}$, but what is relevant here is just the net source strength. For simplicity assume that $[{dS}({E}_{\max },{\boldsymbol{x}}){/{dE}}_{\max }]$ is a power law ${E}_{\max }^{a-1}{{dE}}_{\max }$, (i.e., $S(E\geqslant {E}_{\max })\propto {E}_{\max }^{a}]$) up a to a "maximum" E_max $\equiv \,{E}^{* }$, and vanishes at ${E}_{\max }\gt {E}^{* }$. Ignore, also for simplicity, the possible variation of a with location ${\boldsymbol{x}}$; i.e.,

$$\begin{eqnarray}{dS}({E}_{\max },{\boldsymbol{x}}) & \propto & {E}_{\max }^{a-1}{{dE}}_{\max };\quad {E}_{\max }\leqslant {E}^{* }\\ & \propto & 0;\qquad \qquad \quad {E}_{\max }\gt {E}^{* },\end{eqnarray} \tag{ 2 }$$

then

$$\begin{eqnarray}&&{N}_{p.T}(E,x)\propto {E}^{-2-p}[{{E}^{* }}^{a}-{E}^{a}]/a;\,E\lt {E}^{* }.\end{eqnarray} \tag{ 3 }$$

For secondaries, such as boron, the overall spectrum ${N}_{s,T}$ is the sum of secondaries produced by primaries with the various values of E_max. The contribution from each E_max is proportional to the average grammage $G(E,{E}_{\max })$ traversed by primaries of energy E that are produced with maximum energy E_max. Let us further assume that the relative strength $S({E}_{\max })$ of various points of the Galaxy contributing to N_T is the same for secondaries as it is for primaries; then the total secondary spectrum ${N}_{s,T}$ is given by

$$\begin{eqnarray}{N}_{s,T}(E) & \propto & {\displaystyle \int }_{E}^{\infty }\left[{\displaystyle \int }_{E}^{{E}_{\max }}{N}_{p}(E^{\prime} )K(E^{\prime} ,E)G(E^{\prime} ,{E}_{\max })\right.\\ & & \times [{dS}({E}_{\max }){/{dE}}_{\max }]{dE}^{\prime} \phantom{\rule{}{2.7ex}}]{{dE}}_{\max }.\end{eqnarray} \tag{ 4 }$$

Here $K(E,E^{\prime} )$ is the multiplicity of secondaries at energy E, produced by a collision with primary energy $E^{\prime}$. For a spallation product, if we take the "energy" E to mean energy per nucleon, $K(E,E^{\prime} )$ is, to a good approximation, $\delta (E^{\prime} -E)$, so Equation (4) reduces to

$$\begin{eqnarray}&&{N}_{s,T}\propto {\int }_{E}^{\infty }N(E)G(E,{E}_{\max })[{dS}({E}_{\max }){/{dE}}_{\max }]{{dE}}_{\max }.\end{eqnarray} \tag{ 5 }$$

To simplify further, let us now write G as a function only of E_max but not of E, which allows for the possibility that G is correlated with E_max, while not depending on energy E. Suppose, for example, that $G({E}_{\max })\propto {E}_{\max }^{-\beta }$. Then primary spectra such as carbon would be as before, while boron would be

$$\begin{eqnarray}&&{N}_{s.T}(E,x)\propto {E}^{-2-p}[{{E}^{* }}^{a-\beta }-{E}^{a-\beta }]/[a-\beta ].\end{eqnarray} \tag{ 6 }$$

Now if $a\gt 0$ and $a-\beta \lt 0$, the expression ${E}^{-2-p}[{{E}^{* }}^{a}-{E}^{a}]$ is dominated by the first term in the brackets, because ${E}^{* }\gt E$, while ${E}^{-2-p}[{{E}^{* }}^{a-\beta }-{E}^{a-\beta }]$ is dominated by the second term. If ${E}^{* }\gg E$. then, to a good approximation, Equations (3) and (6) can be written

$$\begin{eqnarray}&&{N}_{p.T}(E,x)\propto {E}^{-2-p},\end{eqnarray} \tag{ 7 }$$

while

$$\begin{eqnarray}&&{N}_{s.T}(E,x)\propto {E}^{-2-p}[{E}^{a-\beta }].\end{eqnarray} \tag{ 8 }$$

That is to say, the primaries can have the same spectrum as the production spectrum, while the secondaries have a steeper one, even though the escape rate is energy independent. This happens when most of the CRs are made in regions where E_max is nearly E*, while most of the secondaries are produced in a minority of regions that just happen to have low E_max and a high value of grammage G.

The case could be made that cool and warm (as opposed to hot) regions, where most of the refractory elements would be locked up in grains and where the density would be highest, contribute a small fraction of primaries and a large fraction of low-energy secondaries. Refractory elements have a higher relative abundance in the cosmic rays than volatile elements (Ellison et al. 1997). This is attributed to the fact that the grains, being charged but massive relative to protons, have a higher rigidity than protons, which guarantees them entry into the Fermi acceleration process. However, in the warm phase, most heavy elements are believed to be locked up in grains. If most of the heavy elements are therefore injected into the diffusive shock acceleration (DSA) process, the refractory elements would be enhanced by far more than the observed factor of ∼4. Thus, it could be argued that most of the cosmic rays do not come from the cool or warm phase, whereas most of the refractories do. (The argument is not airtight, because it could be that the grains are so massive that, even if charged, they are not turned around by the shock, but this would require fine tuning.)

One can apply similar reasoning to secondaries: if CRs are blown out of the disk in collective outflows from multiple supernovae, one might expect that the grammage G they traverse is proportional to the density but inversely proportional to the convection velocity u, while u might anticorrelate as density n as well. Thus, G may depend on density n to a power that exceeds unity—i.e., faster than a linear dependence. This could mean (see below, after further explanation) that most of the secondaries, but only a small fraction of the primaries, are made in the warm phase.

Now consider the quantity E_max. In the following discussion, a specific example is given of how E_max might anticorrelate with G. E_max is limited (e.g., Bell 1978; Drury et al. 1996—hereafter DDK) by the condition that ion-neutral damping does not eliminate the waves that are necessary to confine the accelerated CR near the shock front. I now summarize the calculations of DDK. Suppose that the wave energy density is dissipated by ion neutral damping at the rate of ${\tau }_{d}^{-1}\delta {B}^{2}/4\pi$, where ${\tau }_{d}$ is the wave damping time. When this exceeds the rate of energy gain the waves cannot survive. The rate of wave energy gain per unit volume in the frame of the fluid, ${\dot{U}}_{W}$, due to the force per unit volume exerted by the CR on the incoming fluid, ${ \mathcal F }$, is given by ${v}_{\mathrm{ph}}{ \mathcal F },$ where v_ph is the phase velocity of the waves and ${ \mathcal F }={{dP}}_{\mathrm{CR}}/{dx}$.

Now, at a given energy E, define a quantity ϕ with units of energy flux as

$$\begin{eqnarray}&&\phi =D(E){{dP}}_{\mathrm{CR}}/{dx},\end{eqnarray} \tag{ 9 }$$

where D(E) is the diffusion coefficient of energy E. D(E) is given, in turn, by

$$\begin{eqnarray}&&D(E)\simeq \eta [{B}^{2}/\delta {B}^{2}]{r}_{g}c/3,\end{eqnarray} \tag{ 10 }$$

where r_g is the CR gyroradius, and $\eta \sim 1$ and it is henceforth dropped for convenience. When there is steady state with no escape upstream, $\phi ={{uP}}_{\mathrm{CR}}$ upstream of the shock.

So the condition that the wave growth is positive, $[\delta {B}^{2}/4\pi ]{\tau }_{d}^{-1}\leqslant {v}_{\mathrm{ph}}{{dP}}_{\mathrm{CR}}/{dx}$, can be written as

$$\begin{eqnarray}&&{\tau }_{d}^{-1}\leqslant {v}_{\mathrm{ph}}\phi \left[\displaystyle \frac{4\pi }{{B}^{2}}\right][3{eB}/{Ec}].\end{eqnarray} \tag{ 11 }$$

In order for shock acceleration to work, the diffusive flux ϕ of particles escaping upstream must not exceed the convective flux uP_CR in the shock frame by more than the downstream convective losses ${u}_{+}{P}_{\mathrm{CR},+}$, where ${u}_{+}\equiv {u}_{s}/r$ is the downstream velocity in the shock frame and r is the compression ratio, ${P}_{\mathrm{CR},+}$ is the post shock CR pressure, and, since P_CR at the upstream free escape boundary is probably much less than $({u}_{+}/u){P}_{\mathrm{CR},+}$, we can just say $\phi \leqslant {u}_{+}{P}_{\mathrm{CR},+}$.

$$\begin{eqnarray}&&\phi \leqslant {u}_{+}{P}_{\mathrm{CR},+}\end{eqnarray} \tag{ 12 }$$

whence

$$\begin{eqnarray}&&{\tau }_{d}^{-1}\leqslant 3{v}_{\mathrm{ph}}{u}_{+}{P}_{\mathrm{CR},+}\left[\displaystyle \frac{4\pi }{{B}^{2}}\right]{/r}_{g}c.\end{eqnarray} \tag{ 13 }$$

The damping rate is (DDK)

$$\begin{eqnarray}&&{\tau }_{d}^{-1}=\left(\displaystyle \frac{{\omega }^{2}}{{\omega }^{2}+{({n}_{i}/n{\tau }_{c})}^{2}}\right){\tau }_{c}^{-1}\leqslant {\tau }_{c}^{-1},\end{eqnarray} \tag{ 14 }$$

where

$$\begin{eqnarray}&&{\tau }_{c}^{-1}={\left.8.4\times {10}^{-9}[\displaystyle \frac{T}{{10}^{4}K}\right]}^{0.4}n\,{\mathrm{cm}}^{3}\,{{\rm{s}}}^{-1}\end{eqnarray} \tag{ 15 }$$

(e.g., Parker 1969; Kulsrud & Cesarsky 1971) is the frequency with which a given ion collides with a neutral, and where n is the number density of the neutrals. The time over which the neutrals are dragged along with the ions is $\left(\tfrac{n{\tau }_{c}}{{n}_{i}}\right)$.

Assuming the shock acceleration to be efficient, with an ${E}^{-2}$ differential spectrum, ${P}_{\mathrm{CR}}\sim \rho {u}_{s}^{2}/\mathrm{ln}\left(\tfrac{{E}_{\max }}{{E}_{\min }}\right)$, we use ${P}_{\mathrm{CR}}\,\equiv \rho {u}_{s}^{2}/\mathrm{ln}\lambda$ for brevity, and rewrite Equation (13) as

$$\begin{eqnarray}&&{\tau }_{d}^{-1}\leqslant 3{v}_{\mathrm{ph}}{u}_{+}\left[\displaystyle \frac{4\pi \rho }{{B}^{2}}\right]{u}_{s}^{2}{/r}_{g}c\mathrm{ln}\lambda =\left(\displaystyle \frac{3}{r}\right){v}_{\mathrm{ph}}{u}_{s}^{3}{/r}_{g}{{cv}}_{{\rm{A}}}^{2}\mathrm{ln}\lambda .\end{eqnarray} \tag{ 16 }$$

Note that if the spectrum is steeper than ${E}^{-2}$, $\mathrm{ln}\lambda$ is effectively raised, because the CR pressure at E_max is lowered below $\rho {u}_{s}^{2}/\mathrm{ln}\left(\tfrac{{E}_{\max }}{{E}_{\min }}\right)$.

Assuming the waves are gyroresonant with the CR that generate them, then ${v}_{\mathrm{ph}}={r}_{g}\omega$ where $\omega ={{kv}}_{\mathrm{ph}}$ is the wavefrequency, and, with Equation (14), Equation (16) can be written as

$$\begin{eqnarray}&&\left(\displaystyle \frac{{\omega }^{2}}{{\omega }^{2}+{({n}_{i}/n{\tau }_{c})}^{2}}\right){\tau }_{c}^{-1}\leqslant \left(\displaystyle \frac{3}{r}\right)\omega {u}_{s}^{3}/c{v}_{{\rm{A}}}^{2}\mathrm{ln}\lambda .\end{eqnarray} \tag{ 17 }$$

An upper limit on particle energy, or equivalently, on r_g is obtained only if

$$\begin{eqnarray}&&{\omega }^{2}\gtrsim {({n}_{i}/n{\tau }_{c})}^{2},\end{eqnarray} \tag{ 18 }$$

in which case the condition reduces to

$$\begin{eqnarray}&&{r}_{g}=E/{qB}\leqslant {v}_{\mathrm{ph}}{\tau }_{c}\left(\displaystyle \frac{3}{r}\right){u}_{s}^{3}/c{v}_{{\rm{A}}}^{2}\mathrm{ln}\lambda \equiv {E}_{\max }/{qB}.\end{eqnarray} \tag{ 19 }$$

(If ${\omega }^{2}\ll {({n}_{i}/n{\tau }_{c})}^{2}$, then condition (17) becomes a lower limit on E, which is physically meaningless if particles are unable to reach this lower limit, and in any case empty if they are able to.)

Equation (19) sets an upper limit, E_max, to the energy of

$$\begin{eqnarray}{E}_{{\rm{\max }}} & \sim & {({v}_{{\rm{ph}}}/{v}_{{\rm{A}}})({u}_{{\rm{s}}}/{\rm{c}})}^{3}{(T/{10}^{4}{\rm{K}})}^{-0.4}\\ & & \times {(n/1{{\rm{cm}}}^{-3})}^{-1}{({n}_{i}/1{{\rm{cm}}}^{-3})}^{1/2}(10/{\rm{ln}}{\rm{\Lambda }}){10}^{19}\,{\rm{eV}},\end{eqnarray} \tag{ 20 }$$

which for u = 10³ km s⁻¹ is of order 300 GeV.

So, for example, if CRs are convected out of the Galactic disk at a velocity of ${10}^{-4}{\beta }_{-4}c$, independent of energy (Parker 1969), then they would traverse the disk thickness, ∼100 pc, in about $3{\beta }_{-4}^{-1}$ Myr, and would traverse $\sim 5{\beta }_{-4}(n/1\,{\mathrm{cm}}^{-3})$ g cm⁻². Considering that ${({v}_{\mathrm{ph}}/30\mathrm{km}{{\rm{s}}}^{-1})(B/3\mu {\rm{G}})(T{/10}^{4}{\rm{K}})}^{-0.4}(10/\mathrm{ln}\lambda )$ would be only somewhat less than, if not greater than unity, as its individual factors are of order unity, their spectrum would extend up to at least several hundred ${(n/1{\mathrm{cm}}^{-3})}^{-1}$ GeV or so. So the value of β as defined before Equation (4) is −1.

Now condition (18) itself sets an upper limit on E of

$$\begin{eqnarray}&&E/{qB}\leqslant {v}_{\mathrm{ph}}{\tau }_{c}\left[\displaystyle \frac{n}{{n}_{i}}\right]={v}_{\mathrm{ph}}/\langle \sigma v\rangle {n}_{i},\end{eqnarray} \tag{ 21 }$$

giving

$$\begin{eqnarray}&&{E}_{\max * }\,\\ &&=\,70{[B/3\mu {\rm{G}}]}^{2}{[T{/10}^{4}{\rm{K}}]}^{-0.4}{[{n}_{i}/1{\mathrm{cm}}^{-3}]}^{-3/2}[{v}_{\mathrm{ph}}{/v}_{{\rm{A}}}]\,\mathrm{GeV}.\end{eqnarray} \tag{ 22 }$$

While the phase velocity is usually taken to be the Alfvén velocity v_A, v_ph, in the presence of a strong driving force, the phase velocity may be much greater than v_A (Fiorito et al. 1990), and as high as ${u}_{s}/4\gg {v}_{{\rm{A}}}$ without shutting off shock acceleration. Note that E_max* as defined in Equation (22) is not an upper limit to shock acceleration, but rather an upper limit to the energy range at which shock acceleration would be limited by ion-neutral damping.

To summarize, E_max is either (a) limited by ion-neutral damping to ${E}_{\max }\leqslant {E}_{\max * }=70\left[\tfrac{{v}_{\mathrm{ph}}}{{v}_{{\rm{A}}}}\right]{\left[\tfrac{B}{3\mu {\rm{G}}}\right]}^{2}{\left[\tfrac{T}{{10}^{4}{\rm{K}}}\right]}^{-0.4}{\left[\tfrac{{n}_{i}}{1\mathrm{cm}}\right]}^{-3/2}$, or (b) E_max is not limited at all by ion-neutral damping. It is in any case limited by the size and age of the supernova remnant to about 10⁴ or 10⁵ GeV (Lagage & Cesarsky 1983), which is well below the knee. In the hot phase of the ISM, assuming $T=7\times {10}^{5}$ K and ${n}_{i}=3\times {10}^{-3}$, Equation (22) states that E_max could attain values as high as $5\times {10}^{4}\;{v}_{\mathrm{ph}}{/v}_{{\rm{A}}}$ GeV with ion-neutral damping, allowing eventual escape before severe adiabatic losses set in.

The fact that the spectral index of −2.7 remains nearly constant until well above the established DSA limit claimed by Lagage & Cesarsky (1983) is puzzling, though it has been claimed (J. R. Jokipii 2017, private communication) that shock drift can surpass the DSA limit. But this is in any case a puzzle for any supernova remnant, even if CR trapping is not limited by ion-neutral damping. It may be that, for some unknown reason, the number of supernova with ${E}_{\max }\gtrsim {10}^{4}\,\mathrm{GeV}$ (or the contribution from shock drift) declines as ${E}_{\max }^{a};a\sim -0.7$.

Recall that the total grammage $G\propto (n+{n}_{i})t$ traversed by a cosmic ray before it is convected out of the disk in time t goes as ${n}^{\alpha }$, where, if (a) the material is mostly neutral $(n+\,{n}_{i})$ is non-decreasing with n, and (b) t is non-decreasing with n, then $\alpha \geqslant 1$. Now suppose that dense cloud regions have neutral densities of $n\lesssim {10}^{3}\,{\mathrm{cm}}^{-3}$. Equation (20) then suggests that $G\propto {E}_{\max }^{-\beta }$, where $\beta \geqslant 1$. This implies that even if most of the contribution to the CR we observe is made at ${E}_{\max }\gg 1$ TeV, it is possible that a small fraction is made in dense regions where $n\gtrsim {10}^{2}$ cm⁻³ and E_max ranges from several GeV to over a hundred GeV.

The question of adiabatic losses in incompletely ionized media is pertinent to the relative weights of various contributors to the Galactic CR pool. Although highly ionized media allow acceleration to higher energy, they may, by the same token, require more adiabatic losses once the CR are accelerated. Partially ionized media on the other hand, as E_max decreases with the decreasing expansion velocity of the supernova remnant, release CR at $E\geqslant {E}_{\max }$ that may have been accelerated at earlier stages. This suggests that incompletely ionized regions of the ISM, such as the warm ISM and dense regions of new star formation, may be favorably represented relative to highly ionized regions, such as the hot ISM, simply because they more effectively release the CR produced within them with less adiabatic loss.

Although the hot ISM makes up most of the volume, most of the supernovae may occur in dense regions of newly forming stars, and this is another reason that incompletely ionized parts of the ISM may be favorably represented as CR sources relative to the hot ISM. However, such regions are likely to have a high concentration of young, UV-emitting stars and supernovae, and the interstellar gas in them is likely to make sudden transitions from neutrality to a state of high ionization, so they defy simple parametrization. With that in mind, consider that the density in the star-forming region can be as high as ${10}^{2}\mbox{--}{10}^{3}$ cm⁻³. Suppose a collective blast from multiple supernova forms an expanding superbubble. The Stromgren sphere of photoionization from the young stars is overtaken by the forward shock wave of the expanding superbubble over a timescale of 0.1 Myr (Gupta et al. 2016), and after that, the shock expands into a dense, mostly neutral medium, where, by Equation (20), E_max can be as low as several GeV, and where the grammage traversed by any CRs trapped at the shock can be as high as tens of grams per cm², i.e., a thick target. Because the grammage can be such a steep function of density, steeper than the dependence of $S({E}_{\max })$ on E_max, it is possible that these dense regions spawn most of the secondaries while they contribute little to the primaries.

The question of quantifying the distribution $S({E}_{\max })$, i.e., how the parameters that decide E_max are distributed in the Galactic disk, is beyond the scope of this paper, however, it will be the subject of future research. It is possible that interesting constraints come from the fact that refractory elements have the same slopes as volatiles, which set a limit on how much the warm phase contributes to the Galactic cosmic rays.

The question might arise why the Galactic CR spectral index is −2.7 rather than −2. This is an old question and the usual suggested answers are that most of the volume swept out by supernova shocks occurs when the shocks are no longer at their strongest, and/or that the escape rate from the galaxy is rigidity-dependent. But an additional possibility is that the spectral index a of $S({E}_{\max })$, is negative. This would still beg the question of why the antiprotons have a flatter spectrum than boron. Cowsik & Madziwa-Nussinov (2016) suggest that, at primary energies necessary for antiproton production, the grammage traversed in the nested leaky box (where they hypothesize most of the boron production at lower energies takes place) is less than that in the disk at large. Perhaps the spirit of this suggestion can be adapted to the concept of correlation between E_max and G in view of Equations (20) and (22): these equations suggest that, in the warm phase, beyond sufficiently high energy to make 100 GeV antiprotons, ∼1 TeV, the exponent a, which characterizes the dependence $S({E}_{\max })\propto {E}_{\max }^{a}$, itself changes as a function of E_max, i.e., $S({E}_{\max })$ depends on E_max at ${E}_{\max }\lt {E}_{\max * }$ and that, at higher ${E}_{\max }\gt {E}_{\max * }$, $S({E}_{\max })$ ceases to be E_max-dependent.² If this were the case, and still assuming that escape from the Galaxy is via convection and therefore not significantly energy-dependent, the spectral index of antiprotons would be the same as that of the primaries, as observed. This would predict some flattening of the secondary-to-primary ratio near and above 1 TeV even for spallation secondaries such as boron.

Summary