Gravitino Dark Matter and Flavor Symmetries

In supersymmetric theories without R-parity, the gravitino can play the role of a decaying Dark Matter candidate without the problem of late NLSP decays affecting Big Bang Nucleosynthesis. In this work, we elaborate on recently discussed limits on R-parity violating couplings from decays to antideuterons and discuss the implications for two classes of flavor symmetries: horizontal symmetries, and Minimal Flavor Violation. In most of the parameter space the antideuteron constraints on R-parity violating couplings are stronger than low-energy baryon-number-violating processes. Even in the absence of flavor symmetries, we find strong new limits on couplings involving third-generation fields, and discuss the implications for LHC phenomenology. For TeV scale superpartners, we find that the allowed MFV parameter space is a corner with gravitino masses smaller than O(10) GeV and small $\tan\beta$.


Introduction
In supersymmetric theories, R-parity [1,2] is usually introduced to remove unwanted dimension four operators that would lead to fast proton decay; the renormalizable R-parity violating superpotential is: where the indices are generation indices, i, j, k = 1, . . . , 3, and only antisymmetric combinations of i, j (respectively, j, k) are allowed in λ (respectively, λ ). The first three operators violate lepton number while the last violates baryon number, and both types of operators are involved in proton decay. It is then possible for the proton to be stable if only one type of operators is allowed, leaving B (or L) as an accidental symmetry of the theory [3,4]. This is an aspect of the flavor problems associated with low energy Supersymmetry (SUSY): generic soft terms give large contributions to flavor-changing neutral currents (FCNCs), which can be suppressed by assuming that flavor symmetries govern the structure of the Minimal Supersymmetric Standard Model (MSSM) Lagrangian. Two particularly well motivated types of flavor symmetries are Abelian horizontal symmetries (a la Froggatt-Nielsen [5][6][7]) and Minimal Flavor Violation (MFV) [8,9], according to which the Higgs Yukawa operators are spurions of a SU(3) 5 flavor symmetry under which the full MSSM Lagrangian is invariant.
Under the assumption of these flavor symmetries, definite structures of the RPV couplings are predicted: • with a horizontal U (1) symmetry, the relative structure of the RPV couplings is completely determined by the fermion masses and mixings alone [10][11][12][13]; the baryon number violating (BNV) or lepton number violating (LNV) operators are allowed or forbidden independently. In [10], it was argued that, in order not to disagree with LHC null results, LNV operators should be forbidden altogether when considering sub-TeV SUSY. The BNV couplings λ ijk are written in terms of an overall scale λ 323 and depend on the horizontal charges (we denote the charge of a field by the field symbol itself, Φ ≡ q Φ , and the inter-generational difference between two fields as Φ ij ≡ q Φ i − q Φ j ): • in the MFV framework [8,9], the baryon number violating couplings depend just on tan β and an overall scale factor w , while the lepton number violating operators are naturally suppressed. For tan β 1 we have: (1.5) The coefficient w is not constrained by the flavor structure and should be an O(1) number.
We take these examples as a justification to consider scenarios in which only Baryonic R-parity violation (BRPV) is allowed, while lepton number is conserved (at least to a good approximation). This is the scenario that will be studied in the rest of this paper. It should be noted that in both models (eqs. (1.3) and (1.5)), λ 223 is the largest coupling that does not involve a top in the final state. Implicit in R-parity scenarios is stability of the lightest supersymmetric particle (LSP) which can provide a viable relic Dark Matter (DM) candidate. With a neutralino LSP, this is the usual SUSY WIMP scenario, and problems can arise from late time gravitino decays to the LSP [14], disturbing the predictions of Big Bang Nucleosynthesis (BBN). Alternatively, for a gravitino LSP, it is the NLSP decay to the gravitino that is suppressed by the Planck scale M P and can interfere with BBN. In contrast, R-parity violation allows superpartners to decay directly and quickly into SM particles, 1 solving this problem but at the same time eliminating dark matter candidates from the theory. If, however, the gravitino is the LSP, its decay (see Figure 1) is suppressed by the SUSY breaking scale F (or equivalently, by M P ), 2 by the R-parity violating couplings, and by the superpartner scalẽ m. This naturally allows for lifetimes longer than the age of the universe [16]. Because the 1 As noted above, we consider only baryonic RPV in this paper. Then, late NLSP decays are not a problem for a neutralino or squark NLSPs, but they can be for a stau NLSP. In the second case, heavy superapartners and/or light gravitinos would be needed. 2 For a comprehensive review of gravitino interactions, see Ref. [15]. Figure 1. Gravitino RPV decay; the white vertex marks the 1/M P -suppressed interaction, while the RPV interaction λ ijkũ idjdk is marked by a black dot.
gravitino is unstable, its decays will generate cosmic-rays and high-energy γ-ray emission which can potentially be detected by modern indirect-detection experiments. Given the non-observation of gravitino decay products, we will proceed to set limits on RPV couplings and will compare them to bounds coming from low-energy baryon-number-violating processes (which are especially weak for couplings involving third generation fields). Although this has been studied in the literature, many groups have focused only on the bilinear RPV coupling µ i L i φ u [16][17][18][19][20] with just Refs. [21][22][23][24] discussing the trilinear interactions; weak scale supersymmetry was also frequently assumed. In this paper, we do not set the superpartner scale, we discuss the connection to models with flavor symmetries, which has been unexplored so far, and we show that the limits can be stronger than those from low-energy flavor physics. Following Ref. [25], the decay rate of Figure 1 can be written as in the limit of vanishing masses for the final state particles and at leading order in m 3/2 /m. The lifetime is In this equationm is the common mass scale of the squarks which participate in the process; it is slightly modified in presence of a large hierarchy between different squarks. In particular, the detailed dependence on the squark masses is recovered by substituting the factor 19m 4 3/2 /m 4 in (1.6) with where we have denoted byd j the lightest down squark, and n d is the number of down squarks participating in the process, n = 1, md With the pre-inflationary gravitino abundance washed out during inflation, gravitinos are produced by thermal scattering at reheating and by decays of other fields (such as moduli, or the inflaton). We will show overclosure limits coming from the overproduction of gravitinos, and in the following we will assume T R >m > m 3/2 for the reheating temperature. As a conservative choice, we will assume that the full DM relic abundance is generated in toto at reheating; 3 the second class of processes will just strengthen the overclosure bounds that we are considering. 4 The thermal scattering, with a cross section of order σ ≈ g 2 , overcloses the universe unless (see [15,26] for the precise expression) where the equality holds if the gravitino forms all of the dark matter, as we will assume in the following, and T eq ∼ 1.5 eV is the temperature of matter-radiation equality.
This paper is organized as follows: in Section 2, we review the importance of antideuterons for the indirect detection of dark matter candidates and the coalescence model of antideuteron formation. We then compute and discuss the antideuteron injection spectrum. In Section 3 we derive the upper limits on the RPV coupling λ 223 from the lack of antideuterons and discuss the dependence on the SUSY and SUSY mediation scales. We apply these limits in Section 4, where we discuss the implications for models with flavor symmetries. Finally, we conclude in Section 5.

Antideuterons from Gravitino Decays
Measurements of the cosmic-ray antiproton spectrum by BESS [27][28][29] and PAMELA [30] have provided important constraints on cosmic-ray transport in the galaxy, as well as placed limits on exotic source models such as dark matter annihilations or decays and primordial black hole emission. In the near future, data from the AMS-02 experiment on-board the international space-station will provide the most precise measurements to date. While indirect detection limits on antiprotons currently provide the leading constraints on Rparity violating λ ijk couplings, the production of secondary antiprotons through cosmic-ray spallation processes provides an astrophysical background with considerable uncertainty.
In 2000, Donato et al. [31] proposed new physics searches, specifically neutralino annihilations, using heavier anti-nuclei such as antideuterons, antihelium-3, or antitritium. In contrast to antiprotons, the secondary background for antideuterons is highly suppressed at low energies while gravitino decays produce a peaked spectrum. This happens for three reasons: first, the scattering of cosmic-ray protons with interstellar gas produces (secondary) antiprotons only if the cosmic-ray proton has a total energy above the production threshold E p = 7m p in the galactic rest frame. At these energies, the density of Galactic cosmic-rays is substantial, and the antiproton spectrum below ≈ 5 GeV becomes heavily populated by the astrophysical background. In the case of antideuterons, this threshold is increased to E p = 17m p , where a rapid decrease in the Galactic proton spectrum heavily suppresses the 3 If the universe reheats belowm, gravitinos are not produced thermally. Still, a gravitino relic abundance might be produced by moduli or inflaton decay. 4 On the other hand, these limits can be relaxed with a late entropy injection diluting the relic abundance.
astrophysical background. 5 Second, astrophysical production occurs in a center of mass frame which is highly boosted with respect to the rest of the galaxy, whereas dark matter decays occur at rest. This results in the background spectrum peaking at higher energies than that of dark matter decays, which typically peaks in the non-relativistic regime. Finally, the small binding energy of antideuterons causes them to disintegrate rather than lose energy through inelastic scattering, unlike antiprotons for which such collisions lead to an increased abundance at lower energies. These low-energy antideuteron astrophysical backgrounds are 10-50 times less in magnitude than the primary signals expected from naive thermal dark matter models [32,33], so searches for low energy antideuterons can provide a promising discovery channel for new physics. Heavier elements such as antihelium are even cleaner [34], although the expected signals are too small to be observed with the current generation of experiments as shown by one of the authors (EC) [35] and later independently [36]. For the Baryonic R-parity violating operators under consideration here, constraints on the couplings are currently competitive with those from antiprotons and are expected to improve substantially with the results of AMS-02.
The detection of antinuclei from dark matter decay has been thoroughly investigated in the literature with an emphasis on simple two-body final states such as bb or W + W − (see e.g. [31][32][33][37][38][39][40] for antideuterons and [35,36] for antihelium). Recently, Ref. [24] provided the first antideuteron constraints for gravitinos decaying through a variety of Rparity violating operators. One novel feature of their analysis is the detailed treatment of the Monte Carlo parameters controlling the hadronization model which are tuned to reproduce a wider array of experimental antideuteron production rates. In this paper, much of the same production and propagation framework is used, but we do not vary the hadronization model in order to extract the model-dependent features of gravitino decay, and compare them to standard treatments of decaying dark matter. In doing so, we can provide simple scaling relations which allow BRPV coupling constraints to be easily adapted from future updated measurements and more sophisticated propagation schemes that are presented in the context of two-body decays to heavy quark pairs.
In any process producing antinucleons, it is possible for antiprotons and antineutrons to bind together into a nucleus and produce antideuterons. The traditional formation model, known as the 'coalescence mechanism', was designed to empirically describe nuclei production in heavy-ion collisions based on the phase-space distributions of the constituent nucleons. It possesses a single energy-independent parameter, the coalescence momenta p 0 , and assumes that if any antineutron and antiproton pair have relative invariant 4-momenta 0 , they will fuse and form an antideuteron. The parameter p 0 is then tuned to match collider measurements of d production.
It has long been known that this model cannot accommodate the available data for a single value of the coalescence momenta to better than a factor of ∼ 3. Despite this simplistic model, an improved prescription is largely hindered by limited collider data for production of antideuterons from e + e − collisions at high energies, as well as a lacking 5 The proton spectrum peaks at approximately 10 GeV and subsequently falls off proportionally to understanding of the underlying nuclear formation dynamics. However, recently renewed interest in antideuteron searches have led to at least two important improvements. First, it was pointed out in Ref. [40] that the isotropic nucleon distribution functions used in analytic estimates of formation rates led to an artificial suppression of the d production rate at large center of mass energies. In particular, the jet structure of high-energy showers introduces significant angular correlations between nucleons. One must therefore run Monte Carlo simulations and apply the coalescence mechanism on an event-by-event basis using the simulated phase space distributions of protons and neutrons. Second, it was realized that the antideuteron wave-function is spatially localized to ≈ 2 fm and contributions to the nucleon population from long-lived baryons should be omitted, as they decay at large relative distances from the other particles in the shower. In practice, weakly decaying baryons are then excluded by stabilizing particles with a lifetime τ > 2fm/c with a negligible dependence on this parameter due to the large gap between weak and hadronic timescales.
For our study, we first use Feynrules v2.0 package [41] (using a modified version of the gld-grv [42] and RPV-MSSM [43] model files) to translate our R-parity violating Lagrangian into a UFO format readable by matrix element generators. The matrix elements and phase space for the hard processG →ū idjdk are then generated using MadGraph v5.0 and MadEvent [44]. Finally, these parton level distributions are fed into Pythia 8.1 [45] for showering and hadronization.
In order to fix the coalescence parameter we must choose a value which reproduces a measured rate. As previous studies have noted, a single value of the coalescence momentum cannot simultaneously reproduce rates from different underlying processes such as pp vs e + e − . While this can be slightly improved by tuning the hadronization parameters, we follow previous studies which use electron-positron collisions more likely to resemble a dark matter scenario -i.e. color singlets that are not composite. Following the approach of Refs. [32,33,38], we use e + e − → d measurements from ALEPH at the Z 0 resonance, finding (5.9 ± 1.8 ± .5) × 10 −6 antideuterons per hadronic Z 0 -decay with d momenta 0.62-1.03 GeV/c and polar angle | cos θ| < 0.95 ( [46]). We find a value p A=2 0 = 0.192 ± .030 GeV/c consistent with Refs. [32,33].
In Figure 2 we show the typical antideuteron injection spectra for a gravitino decay of mass m 3/2 =10 GeV, 30 GeV, 100 GeV, 1 TeV, and 10 TeV. In solid lines, we show the spectra from the heaviest accessible channel, which is expected to dominate the decay rate in scenarios with flavor symmetries, while dashed lines show the second heaviest contribution 6 . For comparison, we also show the spectra for a standard dark-matter decay to bb in dotted lines. Shaded bands show the acceptance energies for BESS (red), GAPS (green), the low-energy band of AMS-02 (blue), and the high energy band of AMS-02 (gray). Here we assume that the spectra will be shifted to lower energies as the antideuterons propagate through the heliosphere and shift each band upward in energy due to the Fisk potential φ f = 500 MV acting on a unit electric charge in accordance with the Gleeson & Axford Force Field approximation [47]. The vertical normalization of each energy band is arbitrary and we have slightly offset the BESS band in order to keep the others visible. We note that while the energy range of each experiment is fixed, they are rescaled by a factor m −1 3/2 in these dimensionless coordinates. With the injection spectra now in hand, several observations can be made: 1. Comparing between decay channels, we see that the second lightest quark mass channels have a significantly harder spectrum than the heaviest. For m 3/2 less than a few hundred GeV, these low mass final state channels yield slightly more detectable antideuterons. Such behavior is also evidenced in Ref. [24] where the light quark channels provide the best limits on the trilinear BRPV coupling. Interestingly, this behavior reverses for m 3/2 1 TeV, where the heaviest quark channel dominates by a factor ∼ 20 − 30% over the detectable low energies. One explanation may be the following: Increased jet multiplicity as the 2nd and 3rd generation quarks cascade down to u and d type quarks will divide the gravitino's energy. For low masses, this could sufficiently raise the threshold where heavy channels can consistently form the requisite number of protons and neutrons. When the gravitino mass is very high, each jet will contain energy E m p , and the 3-tiered decay of the top-quark will effectively soften the otherwise harder spectrum.
2. Compared with the 2-body decay to bb, we see a significantly softer spectrum for our gravitino decay in all cases. This results in a mild enhancement in detectable antideuterons of O(50%) for m 3/2 ≈ 50 GeV increasing to a significant factor ≈ 3 above 1 TeV. This is mostly attributed to the higher initial multiplicity of quarks in the final state of the hard process which splits the initial gravitino energy into three final states rather than two. In addition to this, the 3-body phase-space allows the hard jets to occasionally align, and thus increase the probability of a neutron and proton coalescing. In the 2-body case, jets are forced back-to-back for a decay at rest, and are therefore less likely to have cross-jet correlations.
3. The formation model used here is distinct from Ref. [24]. Notably, we use Pythia for hadronization (based on the string fragmentation model) while in Ref. [24], Herwig++ (based on the cluster hadronization model) is used. It has been shown in Ref. [48] that differences between the two different models can lead to substantially different preferred values of the coalescence momentum and variances in the spectrum of anti-deuterons produced. Furthermore, our coalescence momentum is fit to a single data-point at the Z 0 -resonance while the Ref. [24] varies the parameters of the hadronization model in order to reproduce results from e + e − and pp collisions at 50 GeV-7 TeV. We therefore expect to see some level of disagreement at higher gravitino masses. In fact, we do find a significant enhancement in our yield (integrated over the low-energy experimental bands) of around 30% at 50 GeV up to 300% at 1 TeV. As this is an artifact of the underlying hadronization and coalescence model, it occurs independent of the two results enumerated above for which the comparison is based on a common framework.
In order to translate the injection spectra into the observable astrophysical fluxes, one must propagate the d nuclei through two stages: interstellar transport from the position of production to the solar system and modulation of the spectra during through the heliosphere. Interstellar transport of antiprotons and light-nuclei is very well studied but unfortunately still suffers from considerable uncertainties. Our implementation of interstellar propagation follows the standard semi-analytical treatment using the 'two-zone diffusion model' which provides a simplified but good approximation by neglecting energy losses and diffusive reacceleration. The neglect of tertiary processes -i.e. non-annihilating inelastic scatters are treated as annihilations -is well justified and in particular does not redistribute the spectrum of antideuterons toward lower energies. In other words, the Green's function which solves the simplified transport equation is proportional to a δ-function and as a result, the injection spectrum can be factored out of propagation. Similarly, propagation through the solar system in the Force-Field model [47] only introduces an energy dependent scaling and a global shift of the spectrum to 500 MeV lower energy. This implies that event-rates and observable fluxes can be readily compared by computing the ratio of the injection spectra, integrated over the experimental acceptance range. As AMS-02 and GAPS results become available, propagation uncertainties are likely to be reduced based on upper limits on the antiproton spectrum. Recent analyses have already incorporated sophisticated numerical treatments of interstellar and heliospheric propagation and present their findings in terms of annihilation or decay to heavy quarks [32]. Conversion of these rates to case of gravitino dark matter is therefore a simple rescaling according to the integrated ratios of injection read from Figure 2. Our treatment of propagation and conversion of the flux to event rates is completely identical to that presented in Section 3 of Ref. [24] and we therefore omit our own details, pointing the interested reader to the discussion presented there 7 . In the next section we will discuss how the differences between our injection spectra and that of Ref. [24] lead to different limits on the BRPV couplings.

Model-independent limits on RPV couplings
In this section, we will show the model-independent limits on the RPV couplings coming from null observation of antideuterons at the BESS experiment [50] and the future reach of the AMS-02 experiment; as a starting point, we take the limits on λ due to lack of observation of cosmic antideuterons from Ref. [24], where the full gravitino decayG → u idjdk →D + . . . was studied. There, no flavor symmetry was assumed, and the strongest limits were set on the coupling λ 112 (uds channel) in the range 10 GeV ≤ m 3/2 1 TeV, form = 1 TeV.
In Figure 2, we showed the spectra generated by decays of a gravitino in the uds, cbs, and tbs channels for different choices of the gravitino mass: compared to Ref. [24], we find that the number of antideuterons produced in the uds channel in the BESS energy sensitivity range is about 50% higher at low gravitino masses (m 3/2 = 50 GeV), and about a factor of 3 higher at m 3/2 = 800 GeV. The local flux of antideuterons scales as λ 2 . A change in the injection spectra by a factor A therefore strengthens the bounds on λ by a factor of √ A assuming we are in the signal dominated regime (i.e. low relative background flux). We then rescale the 95% confidence level limits on the coupling λ 112 from Figures 7 and 8 of [24] in order to self-consistently employ our injection spectra.
As previously mentioned the number of antideuterons produced in the uds and cbs channels differs only by 20-30% in the energy range accessible by the BESS and AMS-02 experiments. This implies that the limits on the respective couplings differ only by 10-15%. We are particularly interested in the coupling λ 223 , which, according to both eqs. (1.3) and (1.5), is the largest coupling that does not involve a top in the final state. As such, the decay process will proceed through the cbs channel for gravitino masses between the b-quark mass and about 1 TeV. At high gravitino masses (above 1 TeV), the decay involving the coupling λ 323 and a top quark will also be relevant. For such high masses,   Figure 3. 95% CL upper limits on different couplings λ ijk from antideuteron data, compared to the values in Ref. [24]. The region above the lines is excluded. In this graph,m is set at 1 TeV.
the tbs channel gives a slightly higher number of antideuterons when compared to the cbs channel, so that the 95% CL limit on λ 323 will be slightly stronger than the limit on λ 223 at the same scale. As the resulting antideuteron spectrum flattens at the low energies relevant experimentally, we can extend the bounds from Ref. [24] to gravitino masses above 1 TeV and expect no qualitative change in behavior. As we will see in Section 4, because both flavor models predict λ 323 > λ 223 , λ 323 will be the most constrained in this regime. At lower masses, the constraints are strongest for λ 223 . In Figure 3 we compare the limits from Ref. [24] with the ones that will be used in the following. We plot the bounds on the individual couplings λ 112 , λ 323 (λ 223 is degenerate with λ 112 ) as a function of the gravitino mass, with a reference superpartner scalem = 1 TeV.
Provided that the gravitino is the LSP, the injection spectra are independent of the mass of the superpartner involved in the process. Thus, the only dependence on the superpartner scalem is in the hard process that determines the decay rate, as shown in Eq. (1.6). In Figure 4, we show how the limits on λ 223 depend on the superpartner scale: the dot-dashed diagonal lines show the upper limits on λ 223 for given values of m = 1, 10, 100 TeV. The parameter space above each line is ruled out. Alternatively, we fix the ratiom/m 3/2 to approximately 1, 10, 10 2 , 10 3 and show the allowed parameter space. Here the vertical black dashed lines show the upper bounds on the gravitino mass coming from overproduction during reheating. Setting the ratiom/m 3/2 corresponds to setting the SUSY mediation scale M : if the gravitino mass is m 3/2 ≈ F M P and the squark masses arem ∼ F M , we havem/m 3/2 = M P /M . The limits on λ 223 presented in Figure 4 are independent of the flavor structure. As stressed above, for m 3/2 1 TeV, similar limits apply to λ 323 .
It should be noted that the relevant squarks in this process ares R ,c R ,b R ; limits from R-parity conserving LHC searches for first and second generation squarks are above 1 TeV, while forb 1 they are at ≈ 650 GeV [51,52]; without R-parity it is in principle possible for squarks to be significantly lighter than the R-parity conserving constraints. However, we find it a plausible assumption that the superpartners are not hiding at extremely low  masses. Allowing for a little hierarchy betweenb andc,s, we require thatm ≥ 500 GeV for simplicity.

Constraints on models with flavor symmetries
As seen in eqs. (1.3) and (1.5), flavor symmetries constrain the structure of the unknown RPV couplings; limits on one coupling (as λ 223 ) directly translate into limits on all the other couplings (in particular, the largest one, λ 323 ).

Horizontal Symmetries
In models with a horizontal symmetry, as λ 223 is smaller than λ 323 by a factor of 20, the conversion is straightforward: for m 3/2 between 10 GeV and 200 GeV, the cbs channel is predominant in creating antideuterons. Above the top mass, from about 200 GeV to about 1 TeV, the tbs channel contribution grows until it eventually outweighs the cbs channel due to its larger coupling. For m 3/2 1 TeV, the tbs channel gives approximately 20-30% more antideuterons per decaying gravitino than the cbs channel (see Figure 2), with tbs dominating the decays given the larger coupling.
In Figure 5, we present the limits on the largest allowed RPV coupling, λ 323 , which is likely to be the most relevant for LHC phenomenology. On the left, we show limits on λ 323 for given values ofm, while on the right we fix the ratio m 3/2 /m. We also show the future reach of the AMS-02 experiment. An improvement of a factor of 10 is expected across the entire range of gravitino masses.
These limits on RPV couplings can be compared to those found when requiring RPV gives a small contribution to low-energy flavor changing processes, in particular neutronantineutron oscillation or the neutron decay n → Ξ. In [10], one of the authors (AM) showed that the largest RPV coupling is bound to be less than about 10 −2 −10 −3 , depending on the dominant process and the superpartner scale, and independent of the gravitino mass. We see that, apart from m 3/2 30 − 50 GeV, the antideuteron limits from a decaying gravitino DM are stronger than those from low-energy experiments.
Some reference scales should be kept in mind while discussing these limits: • λ 323 = 10 −7 ; in Ref. [53], it was discussed how large R-parity violation would have washed out baryon number in the early universe if the B-violating processes were in equilibrium at a temperature of orderm, and how RPV SUSY at colliders would most likely involve displaced vertices (this was also pointed out in [10] in the context of horizontal symmetries, and in [54]). In Figure 5, the requirement λ 323 < 10 −7 is automatically satisfied for m 3/2 500 GeV for TeV-scale SUSY. For heavier superpartners, or split spectra, it is true for m 3/2 2 − 5 TeV. In other words, for large splittings between the gravitino and the superpartners, the cosmic ray flux from gravitino DM is more constraining than the requirements of having baryogenesis with a large reheating temperature. It should be noted that in baryogenesis scenarios with low reheating temperature [55,56], this bound does not apply as the baryon asymmetry is created after the BNV processes has fallen out of equilibrium (An alternative setting in which baryogenesis is generated by the decay of a meta-stable WIMP was presented in [57].). Given the BICEP2 detection of a tensor-to-scalar ratio r = 0.2 [58], it remains to be seen if such scenarios are still viable.
• λ 323 = 10 −9 : in [10] one of the authors (AM) showed that, in order to evade collider signatures for subTeV SUSY, the lower limit λ 323 > 10 −9 should hold for either a neutralino NLSP (in which case the missing energy signature of R-parity conserving SUSY reappears) or a stop NLSP (for which the the long lived stop hadronizes into Rhadrons and heavy stable charged particle searches would apply). From Figure 5 we can conclude that heavy gravitinos with m 3/2 ∼m > 1 TeV imply λ 323 10 −9 and either give standard R-parity conserving LHC phenomenology or long-lived particles.
• λ 323 = 10 −13 : the lowest scale for which the RPV decay of the NLSP happens before BBN is 10 −13 . We see that this scale is not particularly constrained by Figure 5.
For collider-accessible superpartners, we can conclude that if the coupling λ 323 was measured to be large, it would imply a small gravitino mass.

MFV: A gravitino on the edge
In models with a minimal flavor violating structure [8,9], the only free parameters are the overall scale w , tan β and m 3/2 . The cbs coupling λ 223 is set to λ 223 = 4 × 10 −5 tan 2 βw (4.1) As we will see, we do not need to consider the larger tbs coefficient λ 323 5λ 223 as the limits from the cbs channel are already enough to rule out the large gravitino mass range where the tbs channel (with a top quark in the final state) would dominate.
Setting w ∼ 1, the resulting limits on m 3/2 and tan β are shown in Figure 6. For fixed values ofm, we are forced into a corner with small tan β and/or small m 3/2 . In particular, for LHC-accessible superpartners, the gravitino must be lighter than 50 GeV if tan β ∼ 1, and tan β can be as large as 20 for m 3/2 = 10 GeV. We also note that it is possible to accommodate a 125 GeV Higgs mass in the case of large tan β andm = 10 TeV, as well as in the case ofm = 100 TeV. The AMS-02 experiment will remove a large fraction of this parameter space: form = 1 TeV the limits will be m 3/2 10 GeV, with tan β 5 for m ∼ 10 GeV.
For given values of m 3/2 /m, the only viable options arem = 10 2 m 3/2 andm = 10 3 m 3/2 . In the first case, the gravitino mass should be lower than ∼ 200 GeV for low tan β and below a few tens GeV for larger tan β; in the second case, a larger zone of the parameter space will be explored by AMS-02, but the gravitino is easily overproduced.
If the overall scale factor w was allowed to be 1 a larger region of the parameter space would survive. Given that some couplings are larger than O(10 −7 ), the MFV structure is consistent with high temperature baryogenesis only if w 1. In this case, the limits would scale as √ w and can be relaxed. Still, w cannot be infinitely small, and using the expression (1.5) for the RPV couplings, we avoid the previously discussed limit of λ 323 10 −9 with w 10 −5 . If allowed, a small w should be considered as a tuning of the model.

Conclusions
In this work, we studied how the non-observation of antideuterons cosmic rays places significant constraints on gravitino dark matter in baryonic R-parity violating models with flavor symmetries. We studied a selected number of decay channels and presented limits on the RPV couplings λ 223 and λ 323 , which are almost everywhere stronger than bounds from baryon-number-violating low-energy processes. If flavor symmetries can be used as guides, these are the largest couplings and severe bounds can be cast. While the limits on horizontal flavor symmetries are not as strong, in the minimal flavor violating case the gravitino mass is forced to m 3/2 20 GeV for TeV-scale SUSY. The AMS-02 experiment will be able to reduce this bound below 10 GeV. A suggestive implication (which could hold at least for the MFV scenario) is that the gravitino might be effectively stable, not because of a discrete symmetry such as R-parity, but because decays are not kinematically allowed. 8 Further studies, especially at gravitino masses between 1 and 10 GeV, are needed. In this 8 If the gravitino is lighter than the proton, the proton can decay to it p → K +G . This was considered in Refs. [59,60]), with the most relevant bounds being: range, the best constraints on the RPV coupling will come from antiprotons, positrons and gamma rays. This would also imply a somewhat suppressed mediation scale for SUSY breaking, lower than M P or M GU T , providing a suggestive hint for more new physics at intermediate energies. In a forthcoming publication [61], we are comprehensively exploring all the different decay and detection channels, the uncertainties related to propagation and DM halo profile, as well as the full dependence on the SUSY spectrum.
Note Added While the write-up of this paper was being completed, reference [62] was submitted to the arXiv, which puts similar limits on gravitino DM in the MFV framework by analyzing the antiproton and γ fluxes and has no mention of the SUSY scale. In the present work, the source of the bounds is the lack of observation of antideuterons, which is less sensitive to astrophysical uncertainties. In addition to analyzing other types of flavor symmetries and a larger range of gravitino masses, we extensively discuss sensitivity to the superpartner scale and limits from overclosure.