text:
# Understanding destruction of 𝑛th-order quantum coherence in terms of multi-path interference

## I Introduction

A most profound concept in quantum physics is quantum coherence. The first order version of quantum coherence can be directly manifested by the superposition of two quantum states. This is quite similar to the optical coherence in Young’s double experiment . From the standpoint of the photon- detection theory by Glauber, this simple coherence can be mathematically depicted by the first order correlation function Glauber . However, in quantum mechanics, this first order coherence phenomenon does not sound very marvellous since the same circumstance can also occur in a classical case, such as the optical interference in the above mentioned Young’s double experiment. In fact, by the the first order correlation function only, it is impossible to distinguish the natures of a laser light field and a conventional light field with identical spectral properties. As an effective remedy, Glauber’s $`n`$-th order quantum correlation function, which accounts for various intrinsic ( $`n`$-th order) quantum coherence effects, such as the intensity-intensity correlation measurement in the Hanbury-Brown-Twiss experimentTwiss , was introduced Glauber . Indeed, this function reflects the intrinsically quantum features of coherence beyond the classical analogue.

The quantum coherence embodies the wave nature in the world of microscopic particles. On the other hand, it is very fragile and can easily be destroyed by a “which-path(way)” experiment Rampe ; Umansky , or by an environmental interaction. This phenomenon of destruction of coherence is usually referred to as quantum decoherence Zurek ; Haroche . The phenomenon of losing the coherence described by the first order correlation function is defined as the first order quantum decoherence. Motivated by the considerations in fundamental quantum measurement problem Sun-qm ; Knight , and also by the attempts to preserve quantum coherence of qubits in quantum computing Unruh ; Palma ; Sun-qc , many recent experimental and theoretical investigations have been focused on revealing the physical mechanism of the decoherence problem, e.g., see Sun-adia . According to these studies, this first order decoherence can be roughly understood through the quantum entanglement of the considered system with the environment or the measuring apparatus.

Obviously, this entanglement implies a “which-path(way)” detection Rampe  in the single particle picture. Precisely speaking, in an initial coherent superposition $`|\psi _s=c_n|n`$ , each system state $`|n`$ corresponds to a “path” and many “two path” interferences are reflected in the square norm of the spatial wave function $`x|\psi _s`$.Thus by considering $`x|\psi _s`$ the quantum coherence can be captured to some extent. After the interaction, each “path” is correlated with an environment state $`|e_n`$ to form an entangling state $`|\psi _T=c_n|n|e_n`$. Here, the different states $`|e_n`$ distinguish among the “paths”of different $`|n`$ and thus record the information of each “path”. The interference terms in the spatial intensity $`I(x)=Tr(x|\psi _T\psi _T|x)`$ will disappear when the environment states $`|e_n`$ are completely distinct, i.e., $`e_m|e_n=\delta _{mn}`$. In that case each path is labeled by an environment state .

The above well-known explanation of the first order decoherence in terms of “which-path(way)” detection mechanism is simple but very profound. However, it is not yet clear whether this mechanism can be used to elucidate the $`n`$th order quantum decoherence (n-QDC), the destruction of quantum coherence described by Glauber’s $`n`$th order correlation function. The difficulty is we do not exactly know what are the “paths” and the corresponding “which-path” detection. Most recently, we have touched the second order quantum decoherence problem Zhou-dl  bypassing this difficulty. As a matter of fact, in our treatment, we did not define the concept of “path” directly. The concrete calculation in the ref. Zhou-dl  motivated us to consider the “which-path” picture of the higher order quantum decoherence in general.

In our present investigation, it is crucial to show, we notice that, for a close system, in some cases Glauber’s $`n`$th-order correlation function is the square norm of the $`n`$th-order “multi-particle wave amplitude”, which will be defined later in this paper, while for an open system, it can become a sum of the square norms of the $`n`$th-order “multi-particle wave amplitude” over the states of an environment or an apparatus interacting with this open system. This observation is crucial in our present investigation. As an effective wave function, this multi-particle amplitude can be shown to be a supposition of many generalized “paths” (the multi- particle paths or simply multiple paths). Thus, the higher order coherence may be understood as the interference of multiple particle amplitudes. With this conception a generalized “which-path” detection may be established in terms of “multi-particle paths” as the physical mechanism of higher-order decoherence for some examples.

In section $`2`$, we will briefly explain Glauber’s $`n`$th order quantum coherence in terms of the single and multi photon effective wave functions used in Scully ; Shih . As their generalizations, in section 3 the concepts of multi-particle path and multi-particle wave amplitude are introduced for both close and open bosonic systems. Especially, the $`n`$th order correlation functions of bosonic systems will be studied. In section 4, an interacavity model with two bosonic modes is used to demonstrate $`2`$-QDC as a “which-path” detecting process. In section 5, the exact solution obtained in the appendix is utilized to show the dynamical process of $`2`$-QDC , which is caused by the entanglement with the environment or an apparatus and indeed can be explained as a generalized “which-path” measurement for the explicitly -defined multi-particle paths.

## II The $`n`$th-Order Coherence for Quantized Light Field and Multi-Photon Wave Amplitude

In quantum mechanics, a pure state is a superposition $`|\psi =_kc_k|k`$ of many components $`|k`$ , but the corresponding mixture $`\rho =_k|c_k|^2|kk|`$ can describe the same classical probability distribution $`|c_k|^2.`$ However, $`|\psi `$ and $`\rho `$ represents different quantum realities. Usually it is said that the two components of a pure state is more coherent than those of a mixed one. This coherence property is obviously reflected by the intensity interference of two “paths” corresponding to the two components in the quantum state. For the system of one particle , only single particle property is relevant for this observation. In this sense, only the intensity interference experiment is essential for one particle system. For many particle system, however, there exist many experiments (such as the Hanbury-Brown-Twiss experiment and the intensity-intensity correlation measurement) to show the much richer natures of quantum coherence.

In order to study quantum coherence in many particle system, Glauber introduced the so-called $`n`$-th order quantum correlation function ($`n`$-QCF)

$$G^{(n)}[\alpha _1,t_1;\alpha _2,t_2;\mathrm{};\alpha _n,t_n]=Tr[\widehat{\rho }E_{\alpha _1}^{}(t_1)E_{\alpha _2}^{}(t_2)\mathrm{}E_{\alpha _n}^{}(t_n)E_{\alpha _n}(t_n)\mathrm{}E_{\alpha _2}(t_2)E_{\alpha _1}(t_1)],$$
(1)
for the electro-magnetic field $`E_\alpha (t)`$ in different mode $`\alpha `$. Here, $`E_\alpha (t)`$ is the annihilation operator of mode $`\alpha `$ at time $`t`$ in the Heisenberg picture, $`E_\alpha ^{}(t)`$ is the corresponding conjugate operator, and the density matrix $`\widehat{\rho }`$ represents the initial state of many-mode electro-magnetic field. Of course, this formalism can also be used to study the coherence property for any quantum many-body system.

Furthermore, to describe the higher order coherence, Glauber also defined the $`n`$th-order coherence function in the form footnote

$$g_{\alpha \beta }g[\alpha _1,\alpha _2,\mathrm{},\alpha _n;\beta _1,\beta _2,\mathrm{},\beta _n]=\frac{Tr[\widehat{\rho }E_{\beta _1}^{}(t_1)E_{\beta _2}^{}(t_2)\mathrm{}E_{\beta _n}^{}(t_n)E_{\alpha _n}(t_n)\mathrm{}E_{\alpha _2}(t_2)E_{\alpha _1}(t_1)]}{\sqrt{G[\alpha _1,t_1;\alpha _2,t_2;\mathrm{};\alpha _n,t_n]}\sqrt{G[\beta _1,t_1;\beta _2,t_2;\mathrm{};\beta _n,t_n]}},$$
(2)
Obviously, this coherence function is defined by a ratio of the off-diagonal elements of the reduced multi-time density matrix $`g=(g_{\alpha \beta }:\alpha =(\alpha _1,\alpha _2,\mathrm{},\alpha ),\beta =(\beta _1,\beta _2,\mathrm{},\beta _n))`$to its diagonal ones. It is easily observed from this definition that these off-diagonal elements represent the coherence effect, and each one $`g_{\alpha \beta }`$ correlates two diagonal ones $`g_{\alpha \alpha }`$and $`g_{\beta \beta }.`$ Thus, we can understand the coherence function as measuring the degree of coherence for the two diagonal elements, which correspond to two different “paths”. From this point of view, there is obviously no quantum coherence for a completely mixed state with vanishing off-diagonal elements. In the following discussion, we will give simple examples to illustrate this viewpoint.

Now let us turn to an instance Scully  from which we can see clearly how the higher order quantum coherence effects given by $`2`$-QCF is revealed in the multi-particle picture. For a single photon, the coherent superposition $`|s=\frac{1}{\sqrt{2}}(|1_k+|1_k)`$ of two states with opposite wave vectors $`k`$ and $`k^{}`$ possesses the first order quantum coherence which can be described by the interference fringes

$`G^{(1)}(r,r,t)`$ $`=`$ $`s|E^{}(r,t)E^+(r,t)|s`$ (3)
$`=`$ $`|0|E^+(r,t)|s|^2\mathrm{cos}^2(kr),`$
where

$$E^+(r,t)=E_ka_k\mathrm{exp}(ikri\omega _kt)$$

is the photon field operator with positive frequency for the annihilation operator $`a_k`$ . It should be noticed that the diagonal element $`G^{(1)}(r,r,t)`$ of the first order correlation function is just the square norm of the single photon wave-packet Shih

$`0|E^+(r,t)|s`$ $`=`$ $`{\displaystyle \frac{1}{\sqrt{2}}}\{0|E^+(r,t)|1_k+0|E^+(r,t)|1_k\}`$ (4)
$`=`$ $`{\displaystyle \frac{1}{\sqrt{2}}}(E_ke^{ikri\omega _kt}+E_ke^{ikri\omega _kt})`$
Therefore, $`G^{(1)}(r,r,t)`$ represents the interference between the “two paths” $`0|E^+(r,t)|1_k`$ and $`0|E^+(r,t)|1_k.`$ This just developes the corresponding concept in Young’s double experiment for quantized light field. The above reformulation of the first order quantum coherence implies that “two paths” are necessary for the interference phenomenon.

If we consider the two particle state $`|1_k,1_k^{}`$ with a single component, an interesting situation arises where the first order quantum coherence does not appear, but we can see the second order effect through the second order quantum correlation function

$$G^{(2)}(r_1,r_2,t_1,t_2)=2E_k^4\{1+\mathrm{cos}[(kk^{})(r_1r_2)]\}.$$
(5)
Unlike the case of first order coherence, in this case there do not appear the obvious two or many “paths”. Nevertheless, the interference phenomenon can still be captured in a similar way with the introduction of generalized “path”. We understand the generalized “path” as described by The “two time” correlation function $`G^{(2)}(r_1,r_2,t_1,t_2)`$ $`=|\psi |^2`$ where $`\psi `$ is the “two-photon wave function”

$`\psi `$ $``$ $`\psi (r_1,r_2,t_1,t_2)=00|E^+(r_2,t_2)E^+(r_1,t_1)|1_k,1_k^{}`$ (6)
$`=`$ $`E_k^2e^{2i\omega _kt}[e^{ikr_1+ik^{}r_2}+e^{ik^{}r_1+ikr_2}]`$
$`\psi `$ was invoked as a two photon effective wave-function, and it was also called the biphoton wave packet for the photon field $`E^+(r,t)`$ Shih . Especially, we remark that the biphoton wave packet $`\psi `$ is a coherent superposition of two “probability amplitudes” corresponding to two “two-photon paths”

$$00|E^+(r_2,t_2)|0_k,1_k^{}0_k,1_k^{}|E^+(r_1,t_1)|1_k,1_k^{}$$

and

$$00|E^+(r_2,t_2)|1_k,0_k^{}1_k,0_k^{}|E^+(r_1,t_1)|1_k,1_k^{}$$

Starting from Glauber’s standpoint and proceeding along, we come to the conclusion that a set $`\{G^{(n)}[\alpha _1,t_1;\alpha _2,t_2;\mathrm{};\alpha _n,t_n]|n=1,2,\mathrm{}\}`$of correlation functions, rather than a single one, is indispensable to describe comprehensively the wave-particle dual nature in the quantum world of many particle system. The above example of second order coherence shows that the conception of two-path interference in single photon picture still works with a proper generalization of the concept of “path”. Hence it is quite natural to seek a generalized “which-path(way)” measurement as the mechanism of higher-order quantum decoherence.

## III Multi-Particle Amplitude for free Bosons

In this section the conception of two photon effective wave-function will be generalized. It will be applied to the study of general quantum systems of identical particles. We first discuss the spatially-homogeneous case for the sake of simplicity.

We consider a homogeneous bosonic field with two modes $`|V`$ and $`|H`$. The generalized field operator in time-domain

$`\widehat{\varphi }(t)`$ $`=`$ $`c_V\widehat{b}_Ve^{i\omega _Vt}+c_H\widehat{b}_He^{i\omega _Ht}`$ (7)
$``$ $`c_V(t)\widehat{b}_V+c_H(t)\widehat{b}_H`$
is an annihilation operator with respect to the superposition state

$$|+=c_V^{}|V+c_H^{}|H.$$
(8)
Here , $`b_H`$ and $`b_V`$ are the annihilation operators of the boson system ; $`c_V`$ and $`c_H`$ satisfy the normalization relation $`|c_V|^2+|c_H|^2=1`$. Without loss of generality, we take $`c_V=c_H=1/\sqrt{2}`$ . This means we consider the measurement to detect the polarized boson along the 45o direction in the $`VH`$ plane. We call $`\widehat{\varphi }`$ a “measuring” operator .

Corresponding to $`\widehat{\varphi }`$, the generalized second order correlation function Zhou-dl

$`G^{(2)}`$ $`=`$ $`1_V1_H|\widehat{\varphi }^{}(t_1)\widehat{\varphi }^{}(t_2)\widehat{\varphi }(t_2)\widehat{\varphi }(t_1)|1_V,1_H`$ (9)
$`=`$ $`|0,0|\widehat{\varphi }(t_2)\widehat{\varphi }(t_1)|1_V,1_H|^2|\mathrm{\Psi }(t_1,t_2)|^2.`$
The two time wave function

$$\mathrm{\Psi }(t_1,t_2)=0,0|\widehat{\varphi }(t_2)\widehat{\varphi }(t_1)|1_V,1_H$$

can be understood in terms of the two “paths” picture from the initial state $`|1_V,1_H`$ to the final state $`|0,0`$:

The two “paths” are just associated with the two amplitudes forming a coherent superposition

$`\mathrm{\Psi }(t_1,t_2)`$ $`=`$ $`c_Vc_H\mathrm{exp}(i\omega _Vt_2i\omega _H)+`$ (10)
$`c_Hc_V\mathrm{exp}(i\omega __Ht_2i\omega __Vt_1)`$
Correspondingly, the second order correlation function

$$G^{(2)}=2|c_Vc_H|^2\{1+\mathrm{cos}[(\omega _V\omega _H)(t_2t_1)\}$$
(11)
The above discussion for the second order quantum coherence is applicable to the higher order case. Our arguments in this paper are based on two novel observations: $`a.`$ The generalized field operator $`\widehat{\varphi }=c_n\widehat{b}_n`$ is specified for a quantum measurement about a superposition single particle state $`|\varphi =c_n^{}|n`$. $`b.`$ For a certain initial single component state $`|s_0`$ of $`N`$ particles system, the $`n`$-th order quantum correlation function

$$G^{(n)}(r_1,r_2,\mathrm{},r_n,t_1,t_2,\mathrm{},t_n)=|\psi ^{(n)}|^2$$
(12)
can be written as the norm square of an effective wave function $`\psi ^{(n)}`$, which is just a superposition of many amplitudes.

Let us consider the third order situation as an example. Let the initial state be $`|1_H,2_V`$. Then the generalized third order correlation function

$`G^{(3)}`$ $`=`$ $`2_V1_H|\widehat{\varphi }^{}(t_1)\widehat{\varphi }^{}(t_2)\widehat{\varphi }^{}(t_3)\widehat{\varphi }(t_3)\widehat{\varphi }(t_2)\widehat{\varphi }(t_1)|2_V,1_H`$ (13)
$`=`$ $`|0,0|\widehat{\varphi }(t_3)\widehat{\varphi }(t_2)\widehat{\varphi }(t_1)|2_V,1_H|^2`$
$``$ $`|\mathrm{\Psi }(t_1,t_2,t_3)|^2`$
is a norm square of the two time wave function $`:`$

$`\mathrm{\Psi }(t_1,t_2,t_3)`$ $`=`$ $`0,0|\widehat{\varphi }(t_3)\widehat{\varphi }(t_2)\widehat{\varphi }(t_1)|2_V,1_H`$ (14)
$`=`$ $`\sqrt{2}c_{V}^{}{}_{}{}^{2}c_He^{i\omega _V(t_3+t_2)i\omega _Ht_1}+`$
$`\sqrt{2}c_Hc_{V}^{}{}_{}{}^{2}e^{i\omega _Ht_2i\omega _V(t_3+t_1)}+`$
$`\sqrt{2}c_Hc_{V}^{}{}_{}{}^{2}e^{i\omega _Ht_3i\omega _V(t_2+t_1)}`$
Each term in the above effective wave function is contributed by the corresponding one of the four “paths” from $`|2_V,1_H`$ to $`|0,0`$:         $`|2_V,0_H`$ $`\stackrel{c__V(t_2)}{}`$ $`|1_V,0_H`$     $`\stackrel{c__H(t_1)}{}`$    $`\stackrel{c_V(t_3)}{}`$   $`|2_V,1_H`$ $`\stackrel{c__V(t_1)}{}`$ $`|1_V,1_H`$ $`\stackrel{c_H(t_2)}{}`$ $`|1_V,0_H`$ $`\stackrel{c_V(t_3)}{}`$ $`|0,0`$     $`\stackrel{c_V(t_2)}{}`$  $`\stackrel{c__H(t_3)}{}`$       $`|0_V,1_H`$

In terms of the effective $`3`$-time-wave function defined above, the third order correlation function is explicitly written down:

$`G^{(3)}`$ $`=`$ $`4|c_{V}^{}{}_{}{}^{2}c_H|^2({\displaystyle \frac{3}{2}}+\mathrm{cos}[(\omega _V\omega _H)(t_2t_1)]+`$
$`\mathrm{cos}[(\omega _V\omega _H)(t_3t_1)]+\mathrm{cos}[(\omega _V\omega _H)(t_2t_3)])`$
It shows the quantum interference in the time-domain.

The above analysis is valid only for the case where the considered system is isolated from an environment and not measured by a detecting apparatus-a detector. For our purpose, we need to consider an open system $`S`$ interacting with an environment (reservoir) or a detector $`E`$, and we must extend the concepts of multi-particle(time)-wave functions and the corresponding many-particle paths defined above. To do the generalization , we first invoke the effective field operator

$`\widehat{B}(t)`$ $`=`$ $`U^{}(t,0)\widehat{\varphi }(0)U(t,0)`$ (15)
$`=`$ $`\mathrm{exp}(i\widehat{V}t)\widehat{\varphi }(t)\mathrm{exp}(i\widehat{V}t)`$
where $`\widehat{\varphi }(0)=c_V\widehat{b}_V+c_H\widehat{b}_H`$ has been given by Eq.(7). Then instead of the free time evolution governed by the free Hamiltonian $`H_0`$, we use the evolution operator $`U(t,0)`$ governed by the total Hamiltonian

$$H=H_0+W+H_EH_0+V$$

taking into account the role of the interaction between $`E`$ and $`S`$ . Here, $`H_E`$ is the free Hamiltonian for $`E.`$ If we only consider an ideal quantum decoherence process without dissipation, $`V`$ should possess the nature of quantum non- demolition: $`[H_0,W]=0`$ Sun-qm ; Sun-adia .

Let the states $`|n`$ $``$ $`|n_V,n_H`$ be the common eigen-states of $`H_0`$ and $`W`$ corresponding to the egen-vaslues $`E_n`$ and $`V(n)(n=(n_V,n_H)).`$ If the initial state of the total system is

$$|\psi (0)=|\varphi _s|\varphi _E$$
(16)
where $`|\varphi _s`$ and $`|\varphi _E`$ are some specially-given initial states of $`S`$ and $`E`$ respectively, we can define the effective two-time state vector

$$|\psi _B(t,t^{})=\widehat{B}(t^{})\widehat{B}(t)|\psi (0).$$
(17)
as a reasonable generalization of the effective “two-time wave function” given above. Its norm is just the second order correlation function

$`\psi _B(t,t^{})|\psi _B(t,t^{})`$ $`=`$ $`Tr(\widehat{\rho }(0)\widehat{B}^{}(t)\widehat{B}^{}(t^{})\widehat{B}(t^{})\widehat{B}(t))`$ (18)
$`=`$ $`G^{(2)}[t,t^{},\widehat{\rho }(0)],`$
for the density matrix $`\widehat{\rho }(0)=|\psi (0)\psi (0)|.`$

In fact, due to the non- demolition interaction not resulting in dissipation, the basic dynamic properties of the open system do not change even in the presence of $`E`$ . If we choose $`|\varphi _s=|1_H,1_V,`$ there are only two “paths” from the initial state $`|1_V,1_H`$ to the final state $`|0,0`$ for $`S`$ , and $`|0,0`$ is the unique state which can be reached by the action of $`\widehat{B}(t^{})\widehat{B}(t).`$ Then,

$`\psi _B(t,t^{})|\psi _B(t,t^{})`$ $`=`$ $`{\displaystyle \underset{n,\beta }{}}\psi _B(t,t^{})|n,\beta n,\beta |\psi _B(t,t^{})={\displaystyle \underset{\beta }{}}|0,\beta |\psi _B(t,t^{})|^2`$ (19)
$`=`$ $`{\displaystyle \underset{\beta }{}}|0,\beta |\widehat{B}(t^{})\widehat{B}(t))|1_H,1_V,\varphi _E|^2={\displaystyle }_\beta |\mathrm{\Psi }_\beta (t_1,t_2)|^2`$
Here, the summation ranges over the complete set of states $`|\beta `$ of $`E,`$ and each term in the sum is a norm square of the effective two particle wave function

$$\mathrm{\Psi }_\beta (t_1,t_2)=0,\beta |\widehat{B}(t^{})\widehat{B}(t))|1_H,1_V,\varphi _E$$
(20)
for the open system.

From the above calculations for the second and third order quantum decoherence, we observe that for a specially-given initial state, a higher order correlation function may be explicitly written down as the norm square (or its sum ) of the multi-time-wave function, which is a coherent superposition of several complex components associated with the generalized many-particle paths. It is pointed out that this kind of many-particle path is not a simple-product of single-particle paths, but it can be determined by the specially designed measurement.

## IV Generalized Which-Path Detection in An Intracavity Model

In this section, an intracavity model is presented to demonstrate in multi-particle picture the “which-path” detection associated with higher-order quantum decoherence .

In a recent paper Zhou-dl , we have studied the problem of $`2`$-QDC for a cavity-QED system. The concrete calculation in the ref.Zhou-dl  shows that the $`2`$-QDC effects can indeed be observed in the proposed experiment. But it involves dissipation effect losing energy. However, it is well known that quantum decoherence can still occur for an energy conserving system. So in principle dissipation is not indispensible for the discussion about decoherence effect. For this reason it is natural and interesting to consider pure decoherence process without dissipation. The pure decoherence can be well understood through the quantum entanglement of the considered system with the environment or the measuring apparatus. For a model with pure decoherence process, many concepts (such as multi-particle path and the corresponding which-path detection ) can be made much clearer. Unlike the approximately-solvable model treated in the ref.Zhou-dl , which loses its energy and coherence simultaneously, the model proposed in this section, a bosonic system of two modes interacting with an external system of many harmonic oscillators, does not dissipate its energy. This property makes the model exactly solvable, and as a result the problem of higher-order quantum decoherence can be studied in a straightforward way.

By taking $`\mathrm{}=1,`$ the model Hamiltonian $`\widehat{H}=\widehat{H}_0+\widehat{V}`$ is defined by

$`\widehat{H}_0`$ $`=`$ $`\omega _V\widehat{b}_V^{}\widehat{b}_V,`$ (21)
$`\widehat{V}`$ $`=`$ $`{\displaystyle \underset{j}{}}\omega _j\widehat{a}_j^{}\widehat{a}_j+{\displaystyle \underset{j}{}}[d_V(\omega _j)\widehat{b}_V^{}\widehat{b}_V+d_H(\omega _j)\widehat{b}_H^{}\widehat{b}_H]`$ (22)
$`\widehat{a}_j^{}+\widehat{a}_j),`$
where $`\widehat{H}_0`$ is the free Hamiltonian of the system, $`\widehat{V}`$ the free Hamiltonian $`_j\omega _j\widehat{a}_j^{}\widehat{a}_j`$ of the reservoir (or a detector) plus a non-demolition interaction between the system and the reservoir; and $`\widehat{b}_V^{}(\widehat{b}_V),\widehat{b}_H^{}(\widehat{b}_H)`$ the creation (annihilation) operators for the two modes with frequencies $`\omega _V`$ ($`0`$) and $`\omega _H=0`$. The operators $`\widehat{a}_j^{}(\widehat{a}_j)`$ are creation (annihilation) operators of the reservoir modes of the frequencies $`\omega _j`$. The frequency-dependent constant $`d_H(\omega _j)`$ ($`d_V(\omega _j)`$) measures the coupling constant between $`H`$ ($`V`$) mode and the $`j`$ mode of the reservoir.

The most important feature of the model is the non-demolition condition $`[H_0,V]=0`$. It means that the system does not dissipate energy to the reservoir. On the other hand, the system can leave imprint on the reservoir since, for different number states $`|n_V,n_H`$, there are different interactions

$$\underset{j}{}[n_Vd_V(\omega _j)+n_Hd_H(\omega _j)](\widehat{a}_j^{}+\widehat{a}_j)$$

acting on the oscillator reservoir with different driving forces $`n_Vd_V(\omega _j)+n_Hd_H(\omega _j)`$. When there is only one mode in the external system(reservoir) the whole system can physically be described by an intracavity model: Two mode cavity field interact with a moving wall of the cavity, which is attached to a spring and can be regarded as a harmonic oscillator with a small mass . The fields are coupled to the cavity wall (a moving mirror) by the radiation pressure forces in proportion to the photon numbers $`\widehat{b}_H^{}\widehat{b}_H`$ and $`\widehat{b}_V^{}\widehat{b}_V`$.

For the above introduced model, we can discuss the higher order decoherence problem in the Heisenberg picture by explicitly defining the many-particle “which -path” measurement. The second order coherence is directly determined by the second order correlation function.

$$G[t,t^{},\widehat{\rho }(0)]=Tr(\widehat{\rho }(0)\widehat{B}^{}(t)\widehat{B}^{}(t^{})\widehat{B}(t^{})\widehat{B}(t))$$
(23)
which is defined as a functional of the density operator $`\widehat{\rho }(0)`$ of the whole system at a given time $`0`$. Here, the bosonic field (measuring) operator

$$\widehat{B}(t)=\mathrm{exp}(i\widehat{V}t)(c_H\widehat{b}_H+c_V\widehat{b}_V\mathrm{exp}(i\omega _Vt))\mathrm{exp}(i\widehat{V}t)$$
(24)
is defined for the interacting system. Like the operator defined for the non-interacting system, it also describes a specific destructive quantum measurement  with respect to the polarized states

$$|+=c_H^{}|H+c_V^{}|V$$

where $`c_H`$ and $`c_V`$ satisfy the normalization relation $`|c_H|^2+|c_V|^2=1`$. Without loss of the generality, we take $`c_H=c_V=1/\sqrt{2},`$ considering a specific measuremeant.

To examine whether the macroscopic feature of the reservoir causes the second order decoherence or not , we consider the whole system in an initial state

$$|\psi (0)=|1_H,1_V|\{0_j\},$$
(25)
where $`|\{0_j\}`$ is the vacuum state of the reservoir. Here, we have denoted the general Fock states of the many mode field by $`|\{n_j\}|n_1,n_2,\mathrm{}`$. Because the present discussion concerns the external system interacting the considered system , the conceptions presented in last section must be alternated. Actually, in stead of the effective “two-time wave function”, we use the effective two-time state vector

$$|\psi _B(t,t^{})=\widehat{B}(t^{})\widehat{B}(t)|\psi (0).$$
(26)
Then we can re-write the second order correlation function as

$$G[t,t^{},\widehat{\rho }(0)]=\psi _B(t,t^{})|\psi _B(t,t^{})$$
(27)
It is interesting that the effective state vector can be evaluated as the superposition

$$|\psi _B(t,t^{})=\frac{1}{2}e^{i\widehat{V}(0,0)t^{}}[\mathrm{exp}(i\omega _Vt^{})e^{i\widehat{V}(1,0)t^{}}e^{i\widehat{V}(1,0)t}+\mathrm{exp}(i\omega _Vt)e^{i\widehat{V}(0,0)t^{}}e^{i\widehat{V}(0,1)t^{}}e^{i\widehat{V}(0,1)t}]e^{i\widehat{V}(1,1)t}|\{0_j\}|0_H,0_V$$
(28)
of two components for the two paths from the initial two particle state $`|1_H,1_V`$ to the two particle vacuum $`|0_H,0_V`$ . It should be noticed that the effective actions of the reservoir

$$\widehat{V}(m,n)\underset{j}{}\widehat{V}_j(m,n)=\underset{j}{}\omega _j\widehat{a}_j^{}\widehat{a}_j+\underset{j}{}(d_V(\omega _j)m+d_H(\omega _j)n)(\widehat{a}_j^{}+\widehat{a}_j)$$
(29)
can label the different paths and record the path information in the reservior. Thus , this generalized ”which-path measurement ” leads to the second order quantum decoherence.

The above result clearly demonstrates that, with the presence of the reservoir, the different probability amplitudes ($``$ $`\mathrm{exp}(i\omega _Vt^{})`$ and $`\mathrm{exp}(i\omega _Vt)`$) from $`|1_H,1_V`$ to $`|0_H,0_V`$ entangle with the different states

$$(\frac{1}{2}e^{i\widehat{V}(0,0)t^{}}e^{i\widehat{V}(1,0)t^{}}e^{i\widehat{V}(1,0)t}e^{i\widehat{V}(1,1)t}|\{0_j\}$$

and

$$\frac{1}{2}e^{i\widehat{V}(0,0)t^{}}e^{i\widehat{V}(0,1)t^{}}e^{i\widehat{V}(0,1)t}e^{i\widehat{V}(1,1)t}|\{0_j\})$$

of the reservoir. This is just the physical cause of the second order quantum decoherence. In the following section an explicit calculation of the second order correlation function will be given to illustrate this crucial observation.

## V Dynamic Decoherence in Higher Order Case

After a straightforward calculation the second order correlation function can be expressed in a factorization form :

$$G[t,t^{},\widehat{\rho }(0)]=\frac{1}{2}[1+𝐈m(e^{i\omega _V(tt^{})}\underset{j}{}F_j)]$$
(30)
where each factor

$$F_j=0_j|e^{i\widehat{V}_j(1,1)t}e^{i\widehat{V}_j(0,1)t}e^{i\widehat{V}_j(0,1)t^{}}e^{i\widehat{V}_j(1,0)t^{}}e^{i\widehat{V}_j(1,0)t}e^{i\widehat{V}_j(1,1)t}|0_j0_j|\widehat{u}_j^6(t_6)|0_j$$
(31)
is a two-time transition amplitude of the $`j^{}th`$ mode of the reservoir. Obviously, the term $`_jF_j`$ measures the extent of coherence and decoherence in the second order case. It plays the same role as the decoherence factor of the first order decoherence Sun-qm  . So it is also called the decoherence factor.

In the following, to give the factor $`F_j`$ explicitly, we adopt the Wei-Norman method Wei ; Sun-xiao  to calculate the effective time evolution defined by $`\widehat{u}_j^6(t_6)`$. It can be imagined as an evolution governed by a discrete time-dependent Hamiltonian $`H(t)`$ dominated by $`\widehat{V}_j(1,1),\widehat{V}_j(1,0),\widehat{V}_j(1,0),\widehat{V}_j(0,1),\widehat{V}_j(0,1)`$ and $`\widehat{V}_j(1,1)`$ in six time-intervals $`[t_0=0,t_1=t],[t_1,t_2=2t],[t_2,t_3=2t+t^{}],[t_3,t_4=2t+2t^{}],[t_4,t_5=3t+2t^{}],`$ $`[t_5,t_6=4t+2t^{}]`$ respectively. In the $`k`$-th step of calculation, we take the final state of ($`k1)`$-th step as its initial state. Therefore, we obtain $`\widehat{u}_j^6(t_6)`$ as the sixth step evolution

$$\widehat{u}_j^6(t_6)=e^{g_{1j}^6(t_6)\widehat{a}_{j}^{}{}_{}{}^{}}e^{g_{2j}^6(t_6)\widehat{a}_{j}^{}{}_{}{}^{}\widehat{a}_j}e^{g_{3j}^6(t_6)\widehat{a}_j}e^{g_{4j}^6(t_6)}$$
(32)
Here, $`g_{kj}^6(t_6)(k=1,2,3,4)`$ are the coefficients that can be explicitly obtained. But for the calculation of the $`j`$-th component

$$F_j=\mathrm{exp}[g_{4j}^6(t_6)]$$
(33)
of the decoherence factor, we only need to know $`g_{4j}^6.`$ The detailed discussion in the appendix gives

$`g_{4j}^6(t_6)`$ $`=`$ $`{\displaystyle \frac{2}{\omega _j^2}}[d_H(j)d_V(j)]^2\mathrm{sin}^2[{\displaystyle \frac{1}{2}}\omega _j(t^{}t)]+{\displaystyle \frac{i}{\omega _j^2}}[d_V^2(j)d_H^2(j)]`$ (34)
$`[\omega _j(t^{}t)+2(1\mathrm{cos}(\omega _j[t^{}t])\mathrm{sin}\omega _jt]+(12\mathrm{cos}\omega _jt)\mathrm{sin}(\omega _j[t^{}t])].`$
$`=`$ $`R_j(tt^{})+i\mathrm{\Omega }_j(t,t^{})`$
It is noticed that the real part

$$R_j(tt^{})=\frac{2}{\omega _j^2}[d_H(j)d_V(j)]^2\mathrm{sin}^2[\frac{1}{2}\omega _j(t^{}t)]$$

only depends on the time interval $`t^{}t,`$ but the imaginary part $`\mathrm{\Omega }_j(t,t^{})`$ depends on both $`t^{}`$ and $`t^{}`$ as a two time function.

Because $`R_j`$ can not exceed zero, the norm $`|F_j|=e^{R_j(t)}`$of the factor $`F_j`$ can not exceed one . Then from the arguments about the first order decoherence in our previous works on quantum measurement theory Sun-qm , it is concluded qualitatively that the factorization structure of the decoherence factor

$`F`$ $`=`$ $`{\displaystyle \underset{j=1}{\overset{N}{}}}F_j={\displaystyle \underset{j=1}{\overset{N}{}}}e^{i\mathrm{\Omega }_j(t,t^{})}{\displaystyle \underset{j=1}{\overset{N}{}}}e^{R_j(tt^{})}`$
$`=`$ $`|F|\mathrm{exp}[i\mathrm{\Omega }(t,t^{})]`$
implies the vanishment of the second order correlation in the macroscopic limit that the number $`N`$ of particles making up the reservoir approaches the infinity. This is because

$$|F|=\mathrm{exp}[\underset{j=1}{\overset{N}{}}R_j(t,t^{})]0$$

as $`N\mathrm{}`$ since $`_{j=1}^NR_j(t,t^{})`$ is a diverging series or a monotonously-increasing function of $`tt^{}`$ under some reasonable conditions.

In order to demonstrate the above conclusion quantitatively, we give the numerical results for the second order decoherence for different numbers N of the quantum oscillators. As $`N`$ increase, these results are illustrated in Fig.1. In the numerical calculation, the coupling constants $`\{d_V(j)\}`$ take random values in the domain $`[0.8,1.0]`$, the coupling constants $`\{d_H(j)\}`$ in $`[0.2,0.4]`$, and the frequencies $`\{\omega _j\}`$ in $`[0.5,1.5]`$. The other parameters are given in the caption of the figure.

From both the above graphic illustration given by the numerical calculation and the analytic result given by Eqs.(30,33,34) , we see clearly that the second order correlation function depends not only on the time interval $`t^{}t`$, but also on $`t`$ or $`t^{}`$. In fact it follows from the analytic result

$$𝐈m(e^{i\omega _V(tt^{})}\underset{j}{}F_j)=\frac{1}{2}|F|\mathrm{cos}[\omega _V(tt^{})+\mathrm{\Omega }(t,t^{})]$$
(35)
that the amplitude of the second order correlation function is mainly determined by the time interval $`t^{}t`$, but the phase is determined by the two time parameters . What is more important is that as the number of the quantum oscillators increases, the second order coherence vanishes faster and faster, and the amplitude for quantum revival becomes smaller and smaller. With a reasonable extrapolation, it can be predicted that, when the number of the quantum oscillators approaches the infinity in the macroscopic limit, the second order coherence will vanish in a very short time and thus no quantum revival phenomenon can be observed.

## VI Concluding Remarks

In this paper we first depict the high order quantum coherence of a boson system by introducing the concept of multi-particle wave amplitude. For some cases with the specifically- given initial state , we show that the norm square of multi-particle wave amplitude (or a sum of the norm squares for an open system) gives the high order correlation function. As an effective multi-time wave function, this amplitude can be shown to be a superposition of several “multi-particle paths”. When the environment or an apparatus entangles with them to form a generalized “which-path(way)” measurement, the high order quantum decoherence happens dynamically. Some explicit and general illustrations are presented in this paper to construe our observation. But to prove our conjecture that any high order correlation function indicating the existence of high order quantum coherence can be expressed as the norm square of a properly-defined effective wave amplitude for many particle system, there is still a long way to go. If this conjecture is true under certain general conditions, then there still exist the problem of clarifying these conditions. Moreover, experimental proposals that can be implemented at least in principle are still unavailable.

Our present investigations shed a new light on the understanding of quantum coherence . According to arguments in this paper the high order quantum coherence can also be sculptured as the generalized interference phenomenon by two “multi-particle paths ”, and the intrinsically quantum features (higher order ones) of coherence beyond the classical analogue reflected by the spatial interference of two paths in classical electromagnetic field can be theoretically unified in a framework to embody the wave - particle duality in the quantum world. Therefore, as the essential element, the ”which-path(way)” detection in both the original and the extended versions, naturally provides a complete decoherence mechanism in understanding quantum measurement and the transition from quantum to classical mechanics.

###### Acknowledgements.

This work is supported by the NSF of China and the knowledged Innovation Programme(KIP) of the Chinese Academy of Science. X.F.Liu is acknowledged for his kindly patients in reading and writing this manuscript. The authors also thanks the referee for his some useful suggestions to improve the presentation of this paper.

## Appendixes

In this appendix, the Wei-Norman method Wei ; Sun-xiao  is adopted to calculate the second order decoherence factor $`F_j`$. The calculation is completed in six steps.

During the time period $`[t_{k1},t_k](k=1,2,\mathrm{},6)`$, Let $`W_j^k(t,t^{})`$ be an time evolution dominated by the single particle Hamiltonian

$$\widehat{h}_j^k=\alpha _j^k\widehat{a}_j^{}\widehat{a}_j+\beta _j^k\widehat{a}_j^{}+\gamma _j^k\widehat{a}_j,\{k=1,2,\mathrm{},6\}.$$
(36)
The coefficients $`\alpha _j^k,\beta _j^k,\gamma _j^k`$ and the time intervals $`T_k=t_kt_{k1}`$ take different values in the six different steps:

$`\alpha _j^1=\omega _j,\beta _j^1=\gamma _j^1=d_V(\omega _j)+d_H(\omega _j),T_1=t,`$
$`\alpha _j^2=\omega _j,\beta _j^2=\gamma _j^2=d_V(\omega _j),T_2=t,`$
$`\alpha _j^3=\omega _j,\beta _j^3=\gamma _j^3=d_V(\omega _j),T_3=t^{},`$
$`\alpha _j^4=\omega _j,\beta _j^4=\gamma _j^4=d_H(\omega _j),T_4=t^{},`$
$`\alpha _j^5=\omega _j,\beta _j^5=\gamma _j^5=d_H(\omega _j),T_5=t,`$
$`\alpha _j^6=\omega _j,\beta _j^6=\gamma _j^6=d_V(\omega _j)d_H(\omega _j),T_6=t.`$ (37)
Because of the fact that the four operators $`n_j=\widehat{a}_j^{}\widehat{a}_j`$, $`\widehat{a}_j^{}`$, $`\widehat{a}_j`$ and $`1`$ form a closed algebra - the Heisenberg-Wely algebra, the unitary time evolution operator at each step takes the following form (Wei-Norman theorem)

$$\widehat{u}_j^k(T)=e^{g_{1j}^k(T)\widehat{a}_{j}^{}{}_{}{}^{}}e^{g_{2j}^k(T)\widehat{a}_{j}^{}{}_{}{}^{}\widehat{a}_j}e^{g_{3j}^k(T)\widehat{a}_j}e^{g_{4j}^k(T)}.$$
(38)
for $`T[t_{k1},t_k]`$ in a special sequence. Here coefficients $`g_{sj}^k(T)(s=1,2,3,4)`$ are functions of $`T`$ to be determined. The benefit of the above form is that only the coefficient $`g_{4j}^k(T)`$ is needed in the calculation of the average value at the vacuum state. So we can largely reduce the complexity of our calculation as we need to pay attention only to things concerning $`g_{4j}^k(T)`$ .

Substituting $`\widehat{u}_j^k(T)`$ into the Schr$`\ddot{o}`$dinger equation

$$i\frac{d}{dT}\widehat{u}_j^k=\widehat{h}_j^k\widehat{u}_j^k,$$

we find the coefficients $`g_{sj}^k(T)(s=1,2,3,4)`$ satisfy the following system of equations:

$`{\displaystyle \frac{d}{dT}}g_{2j}^k`$ $`=`$ $`i\alpha _j^k,`$
$`{\displaystyle \frac{d}{dT}}g_{1j}^kg_{1j}^k{\displaystyle \frac{d}{dT}}g_{2j}^k`$ $`=`$ $`i\beta _j^k,`$
$`e^{g_{2j}^k}{\displaystyle \frac{d}{dT}}g_{3j}^k`$ $`=`$ $`i\gamma _j^k`$ (39)
$`{\displaystyle \frac{d}{dT}}g_{4j}^kg_{1j}^ke^{g_{2j}^k}{\displaystyle \frac{d}{dT}}g_{3j}^k`$ $`=`$ $`0`$
Using the results

$$\frac{d}{dT}g_{1j}^k=i\alpha _j^kg_{1j}^ki\beta _j^k,\frac{d}{dT}g_{4j}^k=i\gamma _j^kg_{1j}^k$$
(40)
obtained by simplifying the above system of equations , we get the solution

$`g_{1j}^k(T)`$ $`=`$ $`(g_{1j}^k(t_{k1})+{\displaystyle \frac{\beta _j^k}{\alpha ^k}})e^{i\alpha _j^k(Tt_{k1})}{\displaystyle \frac{\beta _j^k}{\alpha _j^k}},`$
$`g_{4j}^k(T)`$ $`=`$ $`g_{4j}^k(t_{k1})+{\displaystyle \frac{\gamma _j^k}{\alpha _j^k}}(g_{1j}^k(t_{k1})+{\displaystyle \frac{\beta _j^k}{\alpha _j^k}})`$ (41)
$`(e^{i\alpha _j^k(Tt_{k1})}1)+i{\displaystyle \frac{\beta _j^k\gamma _j^k(Tt_{k1})}{\alpha _j^k}}`$
Notice that, to obtain the above result we have used the initial conditions

$`g_{1j}^k(t_{k1})`$ $`=`$ $`g_{1j}^{k1}(t_{k1}),`$ (42)
$`g_{4j}^k(t_{k1})`$ $`=`$ $`g_{4j}^{k1}(t_{k1})`$ (43)
for each step and the initial conditions $`g_{1j}^0(t_0)=g_{4j}^0(t_0)=0`$ for the first step. Then, we obtain a set of iteration equations

$`g_{1j}^k(t_k)`$ $`=`$ $`(g_{1j}^{k1}(t_{k1})+{\displaystyle \frac{\beta _j^k}{\alpha _j^k}})e^{i\alpha ^kT_k}{\displaystyle \frac{\beta _j^k}{\alpha _j^k}},`$
$`g_{4j}^k(t_k)`$ $`=`$ $`g_{4j}^{k1}(t_k1)+{\displaystyle \frac{\gamma _j^k}{\alpha _j^k}}(g_{1j}^{k1}(t_{k1}+{\displaystyle \frac{\beta _j^k}{\alpha _j^k}})`$ (44)
$`(e^{i\alpha _j^kT_k}1)+i{\displaystyle \frac{\beta _j^k\gamma _j^kT_k}{\alpha ^k}}.`$
Iterating six times with different initial conditions and coefficients, the final result of $`g_4^6(t_6)`$ is obtained as Eq.(44).


format: text

text:
# 1 Introduction

## 1 Introduction

Compelling observational evidence suggest that matter in the universe is dominated by “cold dark matter” (CDM), the simplest model of which is one where dark matter particles have interaction strengths at or below the weak scale, and thus interact today only through gravity (collisionless dark matter). The collisionless cold dark matter model has been proven through analytic and numeric simulations to be very successful in explaining the large scale structure of the universe .

However, a series of observations of dark matter structures on the order of $`<1`$ Mpc (e.g., galaxy rotation curves , strong lensing from galaxy clusters , Tully-Fisher relation of spiral galaxies , bar stability in spirals  and deficit of low-mass subhalos in the Local Group ) appear to contradict the prediction of collisionless CDM simulations that dark matter halos form high-density cores, favoring shallower, lower-density dark matter halo cores.

Spergel and Steinhardt  have proposed a variation of CDM in order to alleviate these discrepancies. In their proposed picture, dark matter particles interact with each other strongly and have a large scattering cross section, but negligible annihilation or dissipation.111 If the dark matter particle and antiparticle states were equally populated and annihilations were allowed, these would proceed vigorously, leaving a relic abundance today $`\mathrm{\Omega }_\mathrm{D}10^{37}\mathrm{cm}^2\sigma _{\mathrm{ann}}v/c^1`$ which would be negligible for the case of strong interactions. This implies that the SIDM particles either have strongly suppressed annihilation cross sections, or are asymmetrically populated (as is the case for, e.g., baryons).  Subsequent numerical simulations  have shown that this collisional cold dark matter, or SIDM, does indeed predict halo cores in better agreement with the observations mentioned above, provided

$$0.5<s<6\mathrm{cm}^2\mathrm{g}^1.$$
(1)
where the the parameter $`s=\sigma _{\mathrm{SIDM}}/M_D`$ is the dark matter–dark matter elastic scattering cross section over the dark matter particle mass. Potential difficulties for this model may arise at cluster scales ; one way to address these would be a velocity-dependent cross section, the effect of which we will consider below.

If dark matter particles interact with each other through the strong force, then similar interactions would be expected between them and ordinary baryons as well, with cross sections of the same order . Wandelt et al.  considered this case and found, perhaps surprisingly, that such interactions cannot be excluded by galaxy halo data nor from direct observation by space-borne cosmic-ray detectors, if dark matter particles have a mass larger than $`10^5`$ GeV. Earlier papers  have also considered astrophysical constraints on SIDM both with and without interactions with baryons.

From the particle physics point of view, we choose a model-independent approach. The question we address is “If a massive neutral particle exists, how can we constrain its properties from astrophysical observations?”. We refrain from a detailed discussion of the naturalness of such a particle except to note how quickly model builders proposed several particle physics implementations (e.g. ) of self-interacting dark matter. Some promising candidates, the $`S0`$ and low-mass strangelets, were then ruled out . Nevertheless, there remain intriguing candidates for dark matter which lie in the regions of masses and baryon/self-interaction cross sections we consider, both within supersymmetric extensions of the standard model  or as possible relics from the QCD phase transition within the standard model . Our results place constraints on the energy dependences of the cross sections in these models.

In this paper, we examine the effect of inelastic SIDM–baryon interactions. Specifically, we consider the impact of such interactions on big bang nucleosynthesis (BBN), and on the $`\gamma `$rays produced by interactions of cosmic ray nuclei with dark matter particles. In BBN, the formation of the light elements deuterium, 3He, 4He, and 7Li, depends on the time at which the weak interactions “freeze out,” the free neutron decay rate, and the temperature at which nucleosynthesis occurs, the latter being determined by the photo-dissociation of deuterium. If dark matter can interact strongly with baryons, it is possible that it can destroy light nuclei. Since deuterium has the smallest binding energy and determines the onset of nucleosynthesis, any dark matter deuterium dissociation will change where the onset of nucleosynthesis happens and thus provide a sensitive probe to baryon–dark matter interactions during BBN. We will see that the constraints provided by BBN are quite weak–in other words, for the parameter ranges of interest, SIDM has negligible effect on BBN despite the strength of the SIDM-baryon interaction.

The interaction of cosmic rays with interstellar matter leads to observable signatures in the form of $`\gamma `$rays. Collisions between energetic cosmic ray protons and interstellar hydrogen can produce neutral pions, which in turn decay into two high-energy $`\gamma `$rays; this emission is believed to dominate the observed flux of high-energy photons from the Galactic disk . Including dark matter–baryon interactions allows for similar processes to occur between cosmic ray nuclei and ambient SIDM, thus possibly contributing to the $`\gamma `$ray diffuse background. The $`\gamma `$ray flux from CR-SIDM interactions will be a function of the assumed cross section for the interaction between dark and baryonic matter and can be calculated given the functional form of the Galactic cosmic ray flux and dark matter distribution. The result can then be compared with diffuse $`\gamma `$-ray observations from $`\gamma `$ray telescopes (most recent of which are the EGRET observations, ) and strong upper limits can be placed on the cross sections.

This paper is organized as follows. We will discuss the strong interacting dark matter cross section in §2. In §3 we determine the SIDM effect on BBN and its constraints. In §4 we compute the expected $`\gamma `$-ray flux from cosmic ray - SIDM interaction in the Milky Way and the resulting constraints on the cross section. Other possible constraints are noted in §5. Finally, our findings are summarized and discussed in §6.

## 2 SIDM-Nucleon Interactions

We wish to characterize the inelastic interactions between baryonic nuclei and SIDM particles $`D`$. We will express the cross sections for inelastic interactions in terms of the elastic cross sections. The inelastic behavior is assumed to obey simple scalings with the masses of the reactants  and be energy-independent (but we will use the different constraints below to infer constraints on possible energy dependence of SIDM-baryon interactions).

We begin by considering the SIDM-nucleon interaction, with an energy-independent elastic scattering cross section $`\sigma _{DN}^{\mathrm{elastic}}`$. This need not be the same as the SIDM self-interaction cross section, but was taken to be so in the work of . One can then generalize this to determine the cross section for elastic SIDM-nucleus scattering. Following , we adopt the form

$$\sigma _{DA}^{\mathrm{elastic}}=A^2\left(\frac{\mu (A)}{\mu (p)}\right)^2\sigma _{DN}^{\mathrm{elastic}}$$
(2)
where $`A`$ is the atomic mass number of the interacting nuclide, and $`\mu (A)`$ is the reduced nucleus–dark matter mass. Finally, we will be interested in inelastic scattering, i.e., when internal degrees of freedom are excited in the baryonic component. We will write the (energy-independent) inelastic scattering cross section for channel $`i`$ as

$$\sigma _{DA}^i=\alpha _i\sigma _{DA}^{\mathrm{elastic}}$$
(3)
where $`\alpha _i`$ encodes the relative strength of the inelastic interaction.

A typical strong interaction has a cross section, $`\sigma _{DN}^{\mathrm{elastic}}1`$ barn and it is assumed that a dark matter particle has a mass, $`M_D>1`$ GeV. A useful parameter given these typical values for the cross section and mass is defined as

$$s=\sigma _{DN}^{\mathrm{elastic}}/M_D$$
(4)
usually quoted in cgs units as cm2 g-1. A useful conversion factor for particle physics units is $`1\mathrm{cm}^2\mathrm{g}^1=1.78\mathrm{barn}/(\mathrm{GeV}/\mathrm{c}^2)`$.

## 3 Limiting SIDM with BBN

Primordial nucleosynthesis is the process by which the light elements, consisting mainly of deuterium, 3He, 4He, and 7Li, were produced. Standard BBN is a one parameter theory, depending only on the baryon-to-photon ratio, $`\eta n_\mathrm{B}/n_\gamma =2.74\times 10^8\mathrm{\Omega }_\mathrm{B}h^2`$, where $`\mathrm{\Omega }_\mathrm{B}`$ is the current baryon density relative to the critical density, $`\rho _{\mathrm{crit}}=3H_0^2/8\pi G`$, and $`H_0=100h\mathrm{km}\mathrm{s}^1\mathrm{Mpc}^1`$ is the present value of the Hubble parameter.

Before considering the effect of SIDM, it is useful to recall the basic physics of BBN. At the onset ($`T>1`$ MeV; $`t<1`$ s) of the BBN epoch, the universe was dominated by relativistic particle species, initially consisting of photons, 3 species of light neutrinos, and $`e^\pm `$ pairs; baryons were non-relativistic and essentially only $`n`$ and $`p`$ with no complex nuclei. At this stage, all of these particle species existed in thermal and chemical equilibrium, due to the interaction rates between particles being much faster than the expansion rate of the universe.

At a temperature around an MeV, and time about 1 s, the weak interactions, coupling the neutrinos to the other particle species and keeping the $`n`$ and $`p`$ in equilibrium, become ineffective (interaction rates less than $`H`$). At this time the neutrinos fully decouple from the primordial plasma and the neutron abundance “freezes out”, relative to the proton abundance, with the exception of the occasional free neutron decay. Shortly after, the electrons and positrons start to become non-relativistic, and thus annihilate, producing photons, which in turn heats the plasma relative to the neutrinos. Trace amounts of the light nuclei, deuterium, 3He, 4He, and 7Li exist, with abundances held below nuclear statistical equilibrium levels due to deuterium photo-dissociation by photons in the high-energy tail of the thermal distribution. This “deuterium bottleneck” persists until the universe cools to $`T0.07`$ MeV, when energetic photons are sufficiently rare, at which point deuterium rapidly grows and then is burnt into the other light elements.

We now consider what role CDM plays in these series of events. Depending on the origin of the CDM and its coupling to Standard Model particles, CDM may or may not be in thermal equilibrium with the rest of the cosmic plasma. In the standard CDM scenario, the CDM only interacts weakly with itself and with Standard Model particles, and thus these interactions will freeze out well before the epoch of BBN. Then the dark matter cools rapidly due to Hubble expansion. Throughout BBN, the dark matter is non-relativistic and the Universe is radiation dominated. So in the standard case, dark matter plays no role at all.

What happens in the SIDM case? The dark matter is still non-relativistic, but now elastic collisions between baryons and dark matter will thermalize the dark matter at early times and bring $`T_\mathrm{D}=T_\gamma `$. Dark matter–baryon decoupling will not occur until much later, and we note that this will potentially have interesting effects on structure formation222The tight coupling of the baryon and photon fluids and the coupling between the dark matter and the baryons result in coupling the DM to the photons. This gives pressure to the dark matter fluid and changes the evolution of perturbations until the dark matter decouples from the baryons and then cools adiabatically. We conjecture that this effect leads to a suppression of fluctuation power on perturbation modes which entered the horizon before DM-baryon decoupling. Standard freeze-out arguments show that DM-baryon decoupling occurs at a redshift of $`z2\times 10^4\sigma _{DN}^{\mathrm{elastic}(3/2)}`$, where $`\sigma _{DN}^{\mathrm{elastic}}`$ is measured in barn. The characteristic mass scale of perturbations entering the horizon at this time is of order $`10^9M_{}(\sigma _{DN}^{\mathrm{elastic}})^2`$. Note that this scaling now depends just on the cross section, not on the mass of the particle. This suggests that DM-baryon interactions reduce the dark matter perturbation power on galaxy scales — just those which motivated the study of interacting dark matter in the first place. We plan to study these effects in more detail in a future publication.. However, for our study of BBN elastic scattering is not important except for thermalizing the dark matter. What does affect BBN are inelastic interactions with SIDM particles. Specifically, such collisions can impart internal energy to the nucleus and lead to its dissociation. Given that deuterium has the smallest binding energy and its photo-dissociation determines the onset of nucleosynthesis, we will focus on deuterium breakup; similar processes should occur for other nuclei but are not significant due to the higher binding energies.

We consider the reaction $`D+dD+n+p`$ and neglect the rare, three-body, inverse reaction. In order for this interaction to happen the center-of-mass kinetic energy of the dark matter–deuteron system must be equivalent to at least the binding energy of the deuteron, $`B=2.224`$ MeV. We also assume that the break-up cross section is comparable to the scattering cross section, $`\sigma _{Dd}^{\mathrm{bkp}}=\alpha _1\sigma _{Dd}^{\mathrm{elastic}}`$, where $`\alpha _1`$ represents the efficiency of the break-up reaction and is of order unity. We relate the deuterium–dark matter scattering cross section to the proton–dark matter cross section as defined by eq. (2 ).

In making this interaction applicable for BBN, we need to compute the thermally averaged reaction rate, $`\sigma v`$, where $`v`$ is the relative velocity of the interacting particles. Since for our discussion, we consider a constant cross section with the appropriate energy threshold, this thermal average can be easily found.

$$\sigma v=5.221\times 10^9\alpha _1\sigma _{DN}^{\mathrm{elastic}}\frac{\left(1+\frac{m_p}{M_D}\right)^2}{\left(1+\frac{m_d}{M_D}\right)^{3/2}}\left(1+\frac{B_9}{T_9}\right)T_9^{1/2}\mathrm{exp}\left(\frac{B_9}{T_9}\right)\mathrm{cm}^3\mathrm{s}^1$$
(5)
where $`m_p`$ and $`m_d`$ are the masses of the proton and deuteron respectively, and $`B_9`$ and $`T_9`$ are the binding energy of deuterium and the temperature expressed in units of $`10^9`$ K.333 Throughout this analysis we assume that the dark matter remains non-relativistic during the epoch of BBN, $`M_D>10`$ MeV. If dark matter is relativistic during this time, not only will the reaction rate take on a new form, but now the dark matter will contribute to the expansion of the universe.

The addition of this dark matter reaction adds another channel for deuterium destruction in addition to photo-disintegration. Thus, it will take longer to build up enough deuterium so that nucleosynthesis can occur. To study the effects of the new interaction analytically, we will look at the quasi-static equilibrium (QSE) abundance of deuterium, as laid out by, e.g., , and see how this changes when nucleosynthesis occurs. The abundance of deuterium will be most sensitive to the $`npd\gamma `$, $`d\gamma np`$ and $`dDDnp`$ reactions. The QSE abundance is then given by the ratio of sources to sinks.

$$Y_i^{\mathrm{QSE}}n_i/n_\mathrm{B}=\frac{_{jkl}Y_jY_k[jkil]}{_{jkl}Y_l[iljk]}Y_d=\frac{Y_nY_p[npd\gamma ]}{Y_\gamma [d\gamma np]+Y_D[dDDnp]}$$
(6)
Here $`n_i`$ is the $`i^{\mathrm{th}}`$ particle species number density, $`n_\mathrm{B}`$ is the baryon number density, $`Y_\gamma =1/\eta `$, and $`[jkil]`$ are the reaction rates

$$[jkil]=N_A\rho _\mathrm{B}\sigma _{jkil}v$$
(7)
given by several compilations, where $`N_A=m_u^1`$ is Avogadro’s constant, and $`\rho _\mathrm{B}`$ is the baryon mass density.

The forward and reverse reaction rates are related to each other by detailed balance. This tells us that the ratio of reverse to forward rates is given by the thermal equilibrium distributions of the interacting species. For $`np\gamma d`$ one finds the Saha expression

$$\frac{[d\gamma np]}{[npd\gamma ]}=\frac{Y_n^{\mathrm{EQ}}Y_p^{\mathrm{EQ}}}{Y_d^{\mathrm{EQ}}Y_\gamma ^{\mathrm{EQ}}}=1.40\times 10^5T_9^{3/2}\mathrm{exp}\left(B_9/T_9\right)$$
(8)
Assuming dark matter has a mass $`M_Dm_p`$, the rate is given by the following expression.

$$Y_D[dDDnp]=5.221\times 10^9\alpha _1(\frac{s}{\mathrm{cm}^2\mathrm{g}^1})(\frac{\mathrm{\Omega }_\mathrm{D}}{\mathrm{\Omega }_\mathrm{B}})(\frac{\rho _\mathrm{B}}{\mathrm{g}\mathrm{cm}^3})\left(1+\frac{B_9}{T_9}\right)T_9^{1/2}\mathrm{exp}\left(\frac{B_9}{T_9}\right)\mathrm{s}^1$$
(9)
The ratio of the two deuterium destruction rates is then

$$\frac{Y_D[dDDnp]}{Y_\gamma [d\gamma np]}2.1\times 10^9\eta _{10}\alpha _1(\frac{s}{\mathrm{cm}^2\mathrm{g}^1})\left(\frac{\mathrm{\Omega }_\mathrm{D}}{\mathrm{\Omega }_\mathrm{B}}\right)T_9$$
(10)
assuming $`T_9B_9`$ and where $`\eta _{10}=\eta /10^{10}`$. For reasonable values of $`\alpha _1`$ and $`s`$, we see that this ratio is tiny for the temperatures relevant to BBN. This simply reflects the smallness of $`\eta `$, i.e., the large photon-to-baryon (and photon-to-SIDM) ratio. The smallness of the SIDM contribution will lead to a very weak perturbation to the light elements.

Combining the results of eqs. (8), (9), and (10) and requiring $`\dot{Y}_d=0`$ gives the QSE deuterium abundance. The deuterium bottleneck ends, and the light elements are formed, when this increases to order unity, at the temperature

$$T_9\frac{25.81}{34.18\mathrm{ln}(\eta _{10})\frac{3}{2}\mathrm{ln}(T_9)+2.1\times 10^9\eta _{10}\alpha _1(\frac{s}{\mathrm{cm}^2\mathrm{g}^1})\left(\frac{\mathrm{\Omega }_\mathrm{D}}{\mathrm{\Omega }_\mathrm{B}}\right)T_9}$$
(11)
with $`T_9<1`$.

Thus one can see that with $`\frac{\mathrm{\Omega }_\mathrm{D}}{\mathrm{\Omega }_\mathrm{B}}10`$, one needs an extremely large cross section, $`\alpha _1s>10^7`$ cm2 g-1, to cause a noticeable change in the light element predictions. 4He is the most sensitive to these changes. To derive a limit on SIDM, we will adopt a very generous observational lower bound $`Y_p0.230`$ . Using this, and a full numerical implementation of eq. (9) in the expanded BBN code discussed in , yields a limit $`\alpha _1s1.263\times 10^8\mathrm{cm}^2\mathrm{g}^1`$ for the case $`M_Dm_p`$. Given the rate dependence on dark matter mass and cosmological parameters we can derive an explicit form for the maximum cross section,

$$\alpha _1s<1.263\times 10^9\frac{\left(1+\frac{m_d}{M_\mathrm{D}}\right)^{3/2}}{\left(1+\frac{m_p}{M_D}\right)^2}\left(\frac{\mathrm{\Omega }_\mathrm{B}}{\mathrm{\Omega }_\mathrm{D}}\right)\mathrm{cm}^2\mathrm{g}^1$$
(12)
for a baryon-to-photon ratio $`\eta _{10}=5.8`$ and requiring $`Y_p>0.230`$. The constraints are only slightly different for different assumed $`\eta `$ and helium abundances.

We see, therefore, that BBN cannot place strong constraints on dark matter mass or strong baryon cross section. In other words, energy-independent SIDM is completely compatible with light element constraints and BBN. The weakness of the constraint follows from the smallness of $`\eta `$ and its imposition of the deuterium bottleneck which delays nuclear reactions until temperatures far below the binding energies of the light nuclei. This ensures that SIDM-deuteron interactions will be too weak to allow dissociation, and the lack of SIDM-nucleon bound states ensures that SIDM particles will not otherwise compete with nucleons as the light elements are built up.

Having seen that SIDM-nucleon interactions are compatible with nucleosynthesis in the early universe, we now turn to a present-day consequence of these interactions, namely the production of $`\gamma `$-rays from SIDM collisions with cosmic rays.

## 4 Limiting SIDM with Cosmic Rays

The Galaxy is penetrated by a flux of nonthermal, relativistic nuclei, the cosmic rays. The bulk of these particles (i.e., those with $`E<100`$ TeV) are of Galactic origin, and have kinetic energies typically $`1`$ GeV with a power law energy spectrum

$$\frac{dN}{dE}E^{(\gamma +1)}$$
(13)
where $`\gamma 1.7`$. Inelastic collisions of cosmic ray nuclei (primarily protons) with SIDM particles are detectable in principle by the gamma-ray signature of the decay of the neutral pions produced in such interactions:

$$\begin{array}{ccccc}D+p\hfill & \hfill D+& \mathrm{\Delta }^+\hfill & & \\ & & \hfill & p+\hfill & \pi ^0\hfill \\ & & & & \gamma +\gamma \hfill \end{array}$$

where $`D`$ is the dark matter particle. This reaction is simply the SIDM analog of the usual pion production $`pppp\pi ^0\gamma \gamma `$ reaction, responsible for the bulk of the Galactic contribution to the $`\gamma `$-ray sky .

The Galactic gamma-ray signature for this process is readily calculated. We assume axial symmetry for the cosmic-ray and SIDM spatial distributions, and a relevant dark matter–nucleon cross section $`\sigma _{DN}^{\mathrm{inelastic}}`$ independent of the incoming proton energy. The emissivity $`q_\gamma `$ ($`\gamma `$ray photons per volume per time) is then

$$q_\gamma (r,z)=2n__D(r,z)\sigma _{DN}^{\mathrm{inelastic}}\mathrm{\Phi }_p(r,z).$$
(14)
Here, the factor of 2 is the number of photons per $`\pi ^0`$ decay, $`n__D(r,z)`$ is the number density distribution of the dark matter in the Galaxy, and $`\mathrm{\Phi }_p`$ is the angle-integrated cosmic ray proton intensity for all energies above the pion production threshold,

$$\mathrm{\Phi }_p=4\pi _{E_{\mathrm{min}}}^{\mathrm{}}\varphi _{E,p}𝑑E$$
(15)
with

$$E_{\mathrm{min}}=\left(m_p+m_{\pi ^0}\right)\left[1+\frac{m_{\pi ^0}(2m_p+m_{\pi ^0})}{2M_D(m_p+m_{\pi ^0})}\right]$$
(16)
For $`M_D`$ larger than $`1`$ GeV, the second term in the brackets of eq. (16) is $`1`$ and $`E_{\mathrm{min}}`$ does not change appreciably with $`M_D`$, and we have $`E_{\mathrm{min}}m_p+m_{\pi ^0}`$. Using the parameter $`s=\sigma _{DN}^{\mathrm{elastic}}/M_D`$, and the fact that with typical dark matter halo speeds $`v/c10^31`$, then $`\rho _D=M_Dn_D`$ and we can re-write eq. (14) as

$$q_\gamma (r,z)=2\rho __D(r,z)\alpha _2s\mathrm{\Phi }_p(r,z).$$
(17)
The parameter $`\alpha _2`$ as used in eq. (17), is the product of the ratio $`\sigma _{DpDN\pi }/\sigma _{DN}^{\mathrm{elastic}}`$ and the branching fraction $`f_0`$ for the decay of the $`\mathrm{\Delta }^+`$ to $`p+\pi ^0`$ rather than $`n+\pi ^+`$. From isospin considerations it follows that $`f_0=2/3`$. In the absence of a detailed particle physics model for the interacting dark matter, it is natural to expect that the CR inelastic collision strength is related simply to the breakup collision strength by just this isospin factor, $`\alpha _2=2/3\alpha _1`$, as both reactions are simply exciting internal degrees of freedom in the proton.

The resulting $`\gamma `$ray intensity $`\varphi _{\gamma ,_{DN}}`$ in any given direction will then be the line integral of $`q_\gamma `$ along the line of sight:

$$\varphi _{\gamma ,_{DN}}=\frac{1}{4\pi }_{\mathrm{l}.\mathrm{o}.\mathrm{s}.}q_\gamma (r,z)d\stackrel{}{\mathrm{}}$$
(18)
which, along a line of sight defined by a set of Galactic coordinates ($`\mathrm{},b`$), takes the form

$$\varphi _{\gamma ,_{DN}}=\frac{1}{4\pi }_0^{\mathrm{}}q_\gamma (\sqrt{(R\mathrm{cos}b\mathrm{cos}\mathrm{}L)^2+(R\mathrm{cos}b\mathrm{sin}\mathrm{})^2+(R\mathrm{sin}b)^2},R\mathrm{sin}b)𝑑R$$
(19)
where $`L`$ is the distance of the Sun from the Galactic center, and $`R`$ is the heliocentric distance of a point of Galactocentric distance $`\stackrel{}{r}`$.

A first estimate of the expected $`\gamma `$-ray intensity from cosmic ray - dark matter interactions can be made if we approximate the spatial distribution of the CR flux by a step-function in space: Uniform and equal to the demodulated solar value inside an oblate ellipsoid with a semi-major axis of $`10`$ kpc in the Galactic plane and 1 (5) kpc towards the Galactic poles, and zero outside this ellipsoid. In the limit $`M_D\mathrm{}`$, eqs. (15) and (16) give $`\mathrm{\Phi }_{p,\mathrm{}}=11.8\mathrm{cm}^2\mathrm{s}^1`$, while if we relax our assumption of large $`M_D`$, $`\mathrm{\Phi }_p`$ will scale as $`\mathrm{\Phi }_p=\mathrm{\Phi }_{p,\mathrm{}}\left[1+\frac{m_{\pi ^0}(2m_p+m_{\pi ^0})}{2M_D(m_p+m_{\pi ^0})}\right]^\gamma `$

In addition, we assume that the mass distribution of the dark matter halo can be described as a flattened, non-singular isothermal sphere , $`\rho (r,z)=\rho _cr_h^2/\left(r_h^2+r^2+(z/q)^2\right)`$ , where $`\rho _c`$ is the central density, $`r_h`$ the core radius and $`q`$ the flattening of the halo. In the following calculations we have used $`q=0.8`$, $`\rho _c=0.05\mathrm{M}_{}\mathrm{pc}^3`$ and a value of $`r_h=5\mathrm{kpc}`$, such that the total dark matter mass included in a radius of $`100\mathrm{kpc}`$ is $`10^{12}\mathrm{M}_{}`$. The distance of the Sun from the Galactic center was taken to be $`L8.5\mathrm{kpc}`$.

We have calculated $`\varphi _\gamma `$ in this simple approximation for directions of observation towards the Galactic center $`(\mathrm{}=0^{},b=0^{})`$, towards the anti-center $`(\mathrm{}=180^{},b=0^{})`$ and perpendicular to the Galactic plane $`(b=0^{})`$ and we have compared the results with the corresponding EGRET observational results from Hunter et al  and Sreekumar et al  to derive upper limits for $`\alpha _2s`$. Our results for the limit of large $`M_D`$ are shown in Table 1. The general result for the cross section upper limit will be given by

$$\alpha _2s<\alpha _2s_{\mathrm{}}\left[1+\frac{m_{\pi ^0}(2m_p+m_{\pi ^0})}{2M_D(m_p+m_{\pi ^0})}\right]^\gamma $$
(20)
where $`s_{\mathrm{}}`$ is the value as $`M_D\mathrm{}`$.

Note that for a sufficiently flattened cosmic ray halo (1 kpc minor axis) the strongest constraint comes not from the direction where the $`\gamma `$ray intensity is minimum (the Galactic poles), but from the observations towards the Galactic center, where the $`\gamma `$ray intensity is maximum.

For the derivation of the upper limits presented above, the entire $`\gamma `$ray flux detected by EGRET was attributed to pion decay from cosmic ray–SIDM interactions. This is a good assumption for the direction towards the Galactic poles where the contribution of the $`\gamma `$rays from collisions of cosmic rays and interstellar baryons is not dominant. However, in the Galactic plane most of the observed intensity can be attributed to the presence of interstellar baryonic gas. Thus, the constraints in the directions of the Galactic center and anti-center could be made even more stringent than these conservative limits if one introduces a model for the emission from the interstellar medium.

We note that constraints of this type can be derived for external galaxies as well. In particular, we consider M31, which has the advantages of being nearby, and similar to the Milky Way. M31 has been searched for in EGRET sky maps, but no positive detection has been possible and an upper limit of $`1.6\times 10^8\mathrm{cm}^2\mathrm{s}^1`$ has been placed instead on its diffuse $`\gamma `$ray flux . Given the measured dark mass $`M_D^{\mathrm{M31}}=6\times 10^{10}M_{}`$ in the inner 10 kpc of M31 , we predict a SIDM $`\gamma `$-ray flux of

$$\varphi _\gamma ^{\mathrm{M31}}=2\frac{\mathrm{\Phi }_p^{\mathrm{M31}}\alpha _2sM_D^{\mathrm{M31}}}{4\pi R_{\mathrm{M31}}^2}$$
(21)
where $`R_{\mathrm{M31}}=750`$ kpc is the distance to M31. Following , we estimate the M31 cosmic ray flux $`\mathrm{\Phi }_p^{\mathrm{M31}}`$ by assuming supernovae to be the site of cosmic ray acceleration, and thus taking the cosmic ray flux to be proportional to the supernova rate. This gives $`\mathrm{\Phi }_p^{\mathrm{M31}}=0.45\mathrm{\Phi }_p^{\mathrm{MW}}`$, which allows us to predict $`\mathrm{\Phi }_\gamma ^{\mathrm{M31}}=1.8\times 10^5\alpha _2s\mathrm{cm}^2\mathrm{s}^1`$ (where s is given in $`\mathrm{cm}^2/\mathrm{g}`$). This must be lower than the observed limit, which implies a limit

$$\alpha _2s^{\mathrm{M31}}<9\times 10^4\mathrm{cm}^2\mathrm{g}^1$$
(22)
Thus we see that the M31 $`\gamma `$-ray limits are about an order of magnitude more constraining than those of the Milky Way; however, since they rely on more assumptions regarding the cosmic ray flux in M31, we regard the Milky Way limits as more secure.

Finally, it is very intriguing to note that a recent reanalysis  of the EGRET data finds evidence for a diffuse Galactic $`\gamma `$-ray halo. Specifically, a wavelet analysis was used to identify a $`\gamma `$-ray component at large angular scales, centered on the Galactic center, the intensity of which significantly exceeds the predictions of a model which includes known Galactic sources. One should bear in mind the error budget is large: the inferred halo intensity levels are uncertain to within a factor of $`2`$, and possible systematic effects remain. Also, it is entirely possible that the halo can be explained via inverse Compton scattering of background photons off of energetic cosmic ray electrons .

Nevertheless, it is tantalizing to interpret this diffuse halo in terms of SIDM-baryon interactions. The $`100`$ MeV halo intensity levels toward the poles is quite uncertain but of order $`10^6\mathrm{cm}^2\mathrm{s}^1\mathrm{sr}^1`$, i.e., about 10% of the full EGRET result. In terms of our model, this would amount to a measurement of $`\alpha _2s`$ at about 10% of our limit in Table 1, namely $`\alpha _2s5(1)\times 10^4\mathrm{cm}^2\mathrm{g}^1`$. Future $`\gamma `$-ray observatories such as GLAST can verify the existence of this diffuse Galactic halo, and will be in a position to measure its spectrum. If indeed the flux arises from $`\pi ^0`$ production in SIDM-baryon interactions, we would predict that the spectrum should reflect this origin, in particular, it should be a smooth continuum which shows the “pion bump” feature at $`E_\gamma =m_{\pi ^0}/2=67.5`$ MeV. In addition, similar halo emission should be visible around Local Group galaxies–not only M31, but also the Magellanic Clouds.

## 5 Other Possible Constraints

SIDM-nucleon interactions can have other effects as well. The possibility of DM-baryon interactions contributing to the gamma ray background provokes the question whether some fraction of the $`20`$ unidentified EGRET sources at high Galactic latitudes could be due to dark matter concentrations in the Milky Way halo. A rough estimate shows that typical gamma-ray fluxes of these sources would be compatible with DM clumps of order $`10^3M_{}(R/1\mathrm{kpc})(1\mathrm{cm}^2\mathrm{g}^1/\alpha _2s)`$. While this is an intriguing speculation, the small number statistics of these sources makes distinguishing possible candidates from a background of approximately isotropically distributed AGN a difficult task. A more detailed study will therefore have to await the higher sensitivity and resolution data from future missions such as GLAST.

Another potential (but model-dependent) constraint comes from diffuse radiation. Consider the effect of SIDM particle passage through the interstellar medium of our Galaxy, most of which is hydrogen. Elastic collisions between SIDM particles and the protons occur, with a relative velocity given by $`v200\mathrm{km}/\mathrm{s}10^3c`$, and thus a center-of-mass energy $`E_{\mathrm{CM}}m_pv^2/2=0.2`$ keV. Given that we expect $`m_pM_D`$, it follows that a sizable fraction of this energy will be imparted to the recoiling proton. As a source of heating, this process does not violate observed constraints on the properties of molecular clouds . However, if the scattering process leads to photon emission by the proton, this radiation might be observable. One can show that if each scattering event produces a photon with $`E_{\mathrm{CM}}`$, then the resulting diffuse emission is at or above the level observed in soft X-rays. Of course, it need not be the case that this radiative scattering occurs, but given a particular SIDM model, this constraint should be investigated.

## 6 Discussion and Conclusions

We have considered constraints on strong interactions between dark matter particles and baryons. The big bang nucleosynthesis (BBN) constraints are weak, indicating that interacting dark matter, as it has been recently proposed, is fully consistent with BBN. However, very strong limits come from considering cosmic ray (CR) interactions with the dark matter. For dark matter masses $`M_D>1`$ GeV the ratio $`s\sigma _{DN}/M_D`$ is limited to be less than $`3.9\times 10^3`$cm$`{}_{}{}^{2}/`$g. Constraints from direct detection experiments such as the quantum calorimeter XQC  require $`M_D>10^5`$ GeV. Our limit therefore constrains dark-matter baryon interactions to be less than a hundredth of the self-interaction strength which, simulations suggest, is required to affect the structure of galaxy halos within the self-interacting dark matter (SIDM) scenario. Any particle physics model of SIDM must explain this order of magnitude suppression of the dark matter-baryon interaction.

For the case in which the inelastic cross section is independent of energy, our results are summarized in Figure 1. It is of course possible, and indeed plausible, that SIDM-baryon interactions might be energy-dependent. Thus, we have displayed in Figure 2 the constraints on $`s`$ as a function of the center-of-mass energy at which those constraints apply. Both figures show that energy-independent DM-baryon interactions are ruled out by the cosmic ray constraints, unless the ratio of DM-baryon scattering to DM self-scattering $`\alpha <10^3`$. On the other hand, if the interaction is energy-dependent, Figure 2 constrains the behavior. For example, if we assume that $`sE_{\mathrm{CM}}^n`$, where $`E_{CM}`$ is the center of mass kinetic energy, then we find using our limits from cosmic rays coming from the Galactic center that $`n<0.39`$. Using the M31 constraint steepens this to $`n<0.47`$. In both cases, a $`1/v`$ scaling (or stronger dropoff) is allowed. Thus, our $`\gamma `$-ray constraints provide important information about the energy dependence of a putative SIDM-baryon interaction. If such an energy dependent cross-section were characteristic of dark matter self-scattering as well, SIDM would be less interactive in objects with high velocity dispersion such as clusters of galaxies, possibly softening constraints on SIDM from the statistics of strongly lensed arcs .

We stress that our results do not rely on any detailed assumptions about the underlying particle physics model governing the interactions. In the context of a more detailed particle physics model for SIDM, one could also examine signatures of less generic processes, for example the generation of diffuse emission of soft X-rays from radiative scattering between SIDM and interstellar medium particles.

Finally, we note that elastic SIDM-baryon interactions thermalize the dark matter at early times. This may have potentially important consequences for structure formation. We plan to investigate this issue in more detail in a future publication.

Acknowledgments
We thank Paul Steinhardt for his encouragement and insight, and Jim Buckley for helpful conversations on gamma-ray observations. The work of V.P. was partially supported by a scholarship from the Greek State Scholarship Foundation. The work of B.D.F., V.P., and R.H.C. was supported by National Science Foundation grant AST-0092939.