$\mu$-hybrid Inflation, Gravitino Dark Matter and Stochastic Gravitational Wave Background from Cosmic Strings

We present a successful realization of supersymmetric $\mu$-hybrid inflation model based on a gauged $U(1)_{B-L}$ extension of the minimal supersymmetric standard model, with the soft supersymmetry breaking terms are playing an important role. Successful non-thermal leptogenesis with gravitino dark matter yields a reheat temperature in the range $2 \times 10^{7} \lesssim T_R \lesssim 5 \times 10^{9}$ GeV. This corresponds to the predictions $2 \times 10^{-18} \lesssim r\lesssim 4 \times 10^{-13}$ for the tensor to scalar ratio, and $-2 \times 10^{-6} \lesssim dn_s/d\ln k \lesssim -5 \times 10^{-11}$ for the running of the scalar spectral index. The $B-L$ breaking scale is estimated as $ 6 \times 10^{14}\lesssim M/ \text{GeV}\lesssim 10^{16}$, calculated at the central value of the scalar spectral index, $n_s =0.9655$, reported by Planck 2018. Finally, in a grand unified theory setup the dimensionless string tension parameter associated with the metastable strings is in the range $ 10^{-9} \lesssim G\mu_\text{cs} \lesssim 10^{-6}$ corresponding to a stochastic gravitational wave background lying within the 2$\sigma$ bound of the recent NANOGrav 12.5-yr data.


Introduction
Supersymmetric (SUSY) hybrid inflation [1][2][3][4][5][6][7][8] offers an attractive framework for linking inflation with particle physics models based on grand unified theories (GUTs) [6]. In the minimal models [7,8], the soft supersymmetry breaking terms play an important role in making the predictions for the scalar spectral index n s consistent with the cosmic microwave background (CMB) data [9,10]. Alternatively, the non-minimal terms in the Kähler potential serve a similar purpose [11][12][13][14][15]. With regard to linking inflation with particle physics models, µ-hybrid inflation [16,17] offers an interesting class of SUSY hybrid inflation models where the µ-problem of the minimal supersymmetric standard model (MSSM) is also resolved [16,18]. It is shown in [17] that µ-hybrid inflation based on a renormalizable superpotential and minimal (canonical) Kähler potential leads to split supersymmetry scale with the gravitino mass, m 3/2 5 × 10 7 GeV. However, with non-minimal Kähler potential we can realize µ-hybrid inflation with m 3/2 ∼ 1 − 100 TeV and reheat temperature T R 10 6 GeV [19]. A shifted version of µ-hybrid inflation is investigated in [20] where the monopole problem associated with the breaking of the underlying gauge symmetry can be avoided. A viable model for gravitino dark matter and potentially detectable primordial gravitational waves are among the attractive features of these models. A discussion of successful leptogenesis in µ-hybrid inflation, however, is missing in these papers which is one of the motivations of this paper.
A minimal version of µ-hybrid inflation is considered in the present paper with renormalizable superpotential * E-mail: mansoor@qau.edu.pk and minimal (canonical) Kähler potential. As compared to an earlier treatment of this model in [19], we here allow the soft SUSY breaking mass, M S , to be different from the gravitino mass m 3/2 , as assumed in [8] for standard SUSY hybrid inflation model. This leads to interesting consequences related to the viability of µ-hybrid inflation with gravitino dark matter and successful non-thermal leptogenesis. An adequate range of reheat temperature is obtained while avoiding the gravitino overproduction problem [21,22].
This realization of µ-hybrid inflation is based on a gauged U (1) B−L extension of MSSM, where B and L are the baryon and lepton numbers respectively. The breaking of U (1) B−L gives rise to a topologically stable cosmic string network that is usually constrained from the various experimental bounds. However, here we consider metastable cosmic strings where these bounds can be relaxed. A brief discussion related to the formation of such a metastable string network is presented in a GUT setup based on SO (10). This metastable cosmic string network can decay via the Schwinger production of monopoleantimonopole pairs [23] while generating a stochastic gravitational wave background (SGWB) in a range accessible at the ongoing and future gravitational wave (GW) experiments. We compare our model predictions with the recent bounds from the North American Nanohertz Observatory for Gravitational Waves (NANOGrav) 12.5-yr data [24]. In addition, we highlight the parameter space which is also consistent with gravitino dark matter and successful leptogenesis. For a similar study in no scale inflation see [25].

Supersymmetric µ-hybrid Inflation
The superpotential for hybrid inflation in a U (1) B−L extension of minimal supersymmetric standard model (MSSM) can be written as [26][27][28], where, i, j = 1, 2, 3 are the family indices, κ, λ and λ ij are dimensionless couplings and y are the Yukawa couplings, involving MSSM superfields, The global SUSY minimum occurs at, After the breaking of U (1) B−L gauge symmetry, the last term in eq. (1) gives rise to Majorana mass terms for the right-handed neutrinos, With λ ij (1) and M * ∼ 10 18 GeV, we obtain Majorana masses 10 14 GeV with the gauge symmetry breaking scale M ∼ 10 16 GeV. Therefore, the light neutrino The normalized global SUSY scalar potential V /V 0 as a function of x = |s|/M and y = |φ|/M . masses are naturally generated via the seesaw mechanism. It is also interesting to note that the R−parity which prevents rapid proton decay mediated by the dimension four operators appears as a Z 2 subgroup of U (1) R symmetry. However, the proton is essentially stable due to the global and local symmetries described in table I. A Z 4 subgroup of U (1) R symmetry, consistant with the R charge assignment displayed in table I, is identified as a unique anomaly free discrete symmetry described in [29] that forbids both the µ-term and all dimension four and five baryon and lepton number violating operators in MSSM. The superpotential term relevant for hybrid inflation is [1,2], and the global SUSY F-term scalar potential is given by, where φ, φ, s represents the bosonic components of the superfields Φ, Φ, S respectively. In the D-flat direction, |φ| = |φ|, using eq. (4) and eq. (5), we write the tree level global SUSY potential as, where V 0 = κ 2 M 4 , x = |s|/M and y = |φ|/M . This scalar potential in displayed in fig. 1 where a flat direction (y = 0) with x > 1, suitable for inflation, is clearly visible. As described below, the various important contributions to the scalar potential provide the necessary slope for the realization of inflation in the otherwise flat trajectory. The F-term supergravity (SUGRA) scalar potential is given by, where, φ, s}, z * i being conjugate, and m P 2.4 × 10 18 GeV is the reduced Planck mass. In the present paper we employ the minimal (canonical) Kähler potential given by, and the SUGRA corrections can now be calculated using the above definitions. Along the inflationary trajectory SUSY is broken due to the non-zero vacuum term, V 0 . This generates a mass splitting between the fermionic and the bosonic components of the relevant superfields and leads to radiative corrections in the scalar potential [1]. Another important contribution in the scalar potential arises from the soft SUSY breaking terms [4,5,7]. Including the leading order SUGRA corrections, oneloop radiative corrections and the soft SUSY breaking terms, the scalar potential along the inflationary trajectory (i.e. y = 0) can be written as [7,12,14], where γ ≡ λ/κ and is the one-loop radiative correction function evaluated at the renormalization scale Q, and a is defined as, The last two terms in eq. (9) are the soft SUSY breaking linear and mass-squared terms, respectively, obtained in a gravity-mediated SUSY breaking scheme. The presence of a-term makes the present model a two-field inflation model [30]. However, we assume a suitable initial condition for arg s so that a remains constant during inflation [12]. We further assume the soft mass, M S , of the s field to be different, in general, from the gravitino mass, m 3/2 . As we will see later, this choice provides an extra degree of freedom which yields a relatively wider range of M consistent with the central value of the spectral index n s = 0.966 measured by Planck 2018 [10]. The dimensionless parameter a and the soft mass squared, M 2 S , can have any sign. For standard SUSY hybrid inflation, it is shown in [7] that choosing the negative sign for either soft SUSY breaking term predicts the scalar spectral index n s in good agreement with the central value reported by Planck 2018. We envisage similar results in the present µ-hybrid inflation model.

Inflationary Observables in Slow-roll Approximation
The prediction for the various inflationary parameters are estimated using the standard slow roll parameters defined as, where prime denotes the derivative with respect to x. Note that the extra factor of 1/2 is due to the relation between the canonically normalized real inflaton field, σ ≡ |s|/ √ 2, and the complex field, s. In the slow-roll approximation, the scalar spectral index n s , the tensorto-scalar ratio r and the running of the scalar spectral index α s ≡ dn s /d ln k are given by, The value of the scalar spectral index n s in the ΛCDM model is n s = 0.9665 ± 0.0038 [9]. The amplitude of the scalar power spectrum is given by, where A s (k 0 ) = 2.137 × 10 −9 at the pivot scale k 0 = 0.05 Mpc −1 as measured by Planck 2018 [9]. The relevant number of e-folds, N 0 , before the end of inflation is, where x 0 ≡ x(k 0 ) is the field value at the pivot scale k 0 , and x e is the field value at the end of inflation. As the case may be, the value of x e is fixed either by the breakdown of the slow roll approximation (η(x e ) = −1), or by a 'waterfall' destabilization occurring at the value x e = 1.

Reheating and Leptogenesis
After the end of inflation, the system falls towards the SUSY vacuum and performs damped oscillations about it. The inflaton consists of the two complex scalar fields s and θ = δφ + δφ / √ 2 with the same mass, m inf √ 2κM . The inflaton predominantly decays into a pair of higgsinos ( h u , h d ) and higgses (h u , h d ), each with a decay width, Γ h , given by [31], The other decay mode, via the superpotential couplings leads to a pair of right-handed neutrinos (N ) and sneutrinos ( N ) respectively with equal decay width given as, provided that only the lightest right-handed neutrino with mass M N satisfies the kinematic bound, m inf > 2M N . The relevant Boltzmann equations for the evolution of the total energy density, ρ, of s and θ fields and the radiation energy density, ρ r , are given by, where, With H = 3Γ inf , we define the reheat temperature T R in terms of the inflaton decay width Γ inf , where g * = 228.75 for MSSM. Assuming a standard thermal history, the number of e-folds, N 0 , can be written in terms of the reheat temperature, T R , as [32], Note that the effect of preheating in supersymmetric hybrid inflation is generally expected to be suppressed [33,34]. However, if both the inflaton and waterfall fields are coupled to an additional scalar field, the preheating can be efficient if the inflaton is relatively strongly coupled to this scalar field [33]. In the present case, the electroweak Higgs doublet in the D-flat direction represents such a scalar field and efficient preheating requires λ >> κ. As we only consider λ ∼ κ, the nonperturbative effects via preheating are expected to be suppressed in our case. In addition, preheating associated with fermions is generically expected to be subdominant due to Pauli blocking.
Although sub-dominant, Γ N /Γ h ≤ 1/(3 √ 3γ 2 ) (0.4/γ) 2 < 1, the Γ N channel is important for successful leptogenesis which is partially converted into the observed baryon asymmetry through the sphaleron process [35][36][37]. The washout factor of lepton asymmetry can be suppressed by assuming M N T R . The observed baryon asymmetry is evaluated in term of the lepton asymmetry factor, ε L , where δ ef f is the CP violating phase factor, Γ inf Γ h and, assuming hierarchical neutrino masses, ε L is given by [38], Here, the atmospheric neutrino mass squared difference is ∆m 2 31 ≈ 2.6 × 10 −3 eV 2 and H u 174 GeV in the large tan β limit. For the observed baryon-to-photon ratio, n B /n γ = (6.12 ± 0.04) × 10 −10 [39], the bound on |δ ef f | ≤ 1 along with the kinematic bound, m inf ≥ 2M N , translates into a bound on the reheat temperature, for γ ≥ 1. In order to attain the minimum possible reheat temperature we set γ = 1. This bound from successful leptogenesis, i.e T R 2 × 10 7 GeV, is represented by the gray dashed line in figs. 2-4. An underproduction of leptogenesis is assumed for T R < 2 × 10 7 GeV with a reduction in M N .

Cosmic String Constraints
Cosmic Strings (CSs) are produced at the end of inflation with implications related to anisotropies in the CMB and the production of stochastic gravitational waves (SGWs). The predictions related to SGWs are discussed in section 9. The strength of the string's gravitational interaction is expressed in terms of the dimensionless string tension, Gµ cs , where G = 1/8πm 2 P and µ cs is mass per unit length of the string. The CMB bound on the CS tension is [9,10], The quantity µ cs , can be written in terms of the U (1) B−L gauge symmetry breaking scale M [40], where g = 0.7 for MSSM. Requiring M 10 16 GeV, Gµ cs 10 −6 and this is possible with a metastable CS network as described in [23]. This possibility not only circumvents the CMB bound on Gµ cs , it can also evade other bounds coming from LIGO O3 [41]. A possible realization of a metastable CS network in a GUT setup based on SO(10) model is described in section 9.

Results and Discussion
The seven parameters of the present model, κ, M, am 3/2 , M 2 S , x 0 , x e , and M N , are constrained by five conditions, namely, the amplitude of scalar power spectrum, A s (k 0 ) = 2.137 × 10 −9 eq. (14), the scalar spectral index, n s = 0.9665, the end of inflation by the waterfall, x e = 1, the number of e-folds, N 0 , defined in eq. (15) and given in terms of T R by eq. (21), and finally the observed value of baryon-to-photon ratio, n B /n γ = 6.12 × 10 −10 eq. (22). This leaves two independent parameters to freely vary which can be taken to be am 3/2 and M 2 S . Keeping one parameter fixed, the second can be varied as depicted in figs. 2a and 2b. In fig. 2a, we vary am 3/2 for various values of M 2 S in the range 2 × 10 12 GeV 2 to −2 × 10 17 GeV 2 . In fig. 2a, the curve with M S ∼ 0 has already been discussed in [19] with a minimal (canonical) Kähler potential. On the other hand, the region with |M 2 S | = 0 is the new parametric space explored in this paper. However, see [8] where this region is explored in the standard hybrid inflation model (γ 1). Similarly, in fig. 2b we vary M 2 S for fixed values of m 3/2 lying in the range from 0 to 730 TeV (0 to 155 GeV) for a = 1 (−1).
In accordance with the outcome of SUSY hybrid inflation model [7,8], at least one of the two parameters, M 2 S or am 3/2 , is expected to be negative in order to realize the red-tilted scalar spectral index consistent with Planck-2018 data. In the present model, the scalar spectral index with x 0 ∼ 1 can be written as, In the limit where M 2 S term is dominant in the above expression we obtain,  fig. 2b with fixed values of am 3/2 > 0, it is noted that the soft SUSY breaking and radiative correction terms are comparable in (x 0 ) in order to satisfy the constraint of A s given in eq. (14). This leads to M ∝ κ −3 behavior of curves in the lower region of fig. 2b.
The four boundary curves in figs. 2-4 are respectively described by Gµ cs = 10 −6 , 10 −11 , T R = 10 10 GeV and x 0 = 1.0001. We do not consider larger values of reheat temperature, T R 10 10 GeV, which are usually constrained by the gravitino overproduction problem and allow up to 0.01% difference between x 0 and x e = 1, since for smaller values of κ the corresponding field value behavior obtained from eqs. (16) and (20). Further, the curves with fixed values of inflaton mass, m inf √ 2κM , ranging from 2 × 10 9 GeV to 8 × 10 11 GeV are shown in fig. 3b and are consistent with M ∝ κ −1 behavior.
The predicted range of the tensor to scalar ratio with tiny values, r ∼ 5 × 10 −11 − 5 × 10 −21 , is shown in fig. 4a where the various curves with constant values of r follow M ∝ κ −1/2 behavior as can be obtained from eq. (14), Finally, the relevant expression of α s in the slow-roll approximation is given by, The predicted range −α s ∼ 2.7 × 10 −3 − 1.2 × 10 −13 is shown in fig. 4b where the curves with constant values of α s follow M ∝ κ, that can be obtained from eq. (30) assuming a dominant contribution from the radiative correction with F (x 0 ) ∝ κ −2 . The predicted ranges of r and α s with tiny values are consistent with the underlying assumption of the ΛCDM model.

Gravitino Dark Matter
Following [17,19,20,42], an interesting realization of a stable gravitino as a cold dark matter (DM) candidate is presented here. The relic abundance for stable gravitinos is described in terms of the reheat temperature T R and the gluino mass, mg, as [43][44][45][46], This expression contains only the dominant QCD contributions for the gravitino production rate. The electroweak contribution [44] is expected to be relatively suppressed for our analysis. The observed DM abundance requires, Ω 3/2 h 2 ∼ 0.12 [9], which allows us to write the gravitino mass in terms of the reheat temperature for a given value of gluino mass. The gray dotdashed curves in figs. 2-4 are derived from eq. (31) by taking into account mg m 3/2 and the LHC bound on the gluino mass, mg > 2.2 TeV [47]. The region to the left of these curves describes the gravitino DM in totality. It covers the values of gravitino mass in the range, m 3/2 ∼ 6 GeV − 63 TeV with reheat temperature up to 6 × 10 9 GeV. Hence, the viable parameter space compatible with both DM and leptogenesis lies in the region bounded by the gray dashed and dot-dashed curves of figs. 2-4. Assuming an underproduction of leptogenesis, the region left of the gray dot-dashed curve (see figs. [2][3][4] is also compatible with gravitino DM. With LSP gravitino the next to lightest supersymmetric particle (NLSP) X can decay into SM particles and gravitino. In this case the lifetime of X should be less than 1 sec in order to keep the successful predictions of big bang nucleosynthesis (BBN) intact. The decay rate for X −→ ψ µ γ is given by [48], where θ W is the electroweak mixing angle. For the NLSP lifetime, τX 1 sec, eq. (32) yields a lower bound on For mX ∼ 1 TeV (10 TeV), we obtain 1 GeV m 3/2 < 1 TeV (10 GeV m 3/2 < 10 TeV). This again confirms that the region bounded by the gray dot-dashed and dashed curves is consistent with a gravitino DM and successful leptogenesis.
In order to suppress the washout effects in non-thermal leptogenesis we usually require M N to be somewhat larger than T R . This can be achieved in the present model by exploiting the freedom in the allowed range of the CP-violating phase factor, δ ef f ≤ 1. We can choose any value of the lightest right-handed neutrino mass, M N , lying in the range, T R ≤ M N ≤ m inf /2. For example with M N = m inf /4 we obtain 10 M N /T R 200 consistent with gravitino DM scenario.
For the unstable long-lived gravitino, m 3/2 < 25 TeV, we have to take into account the BBN bounds on the reheat temperature [48,49], which yields 5 TeV m 3/2 25 TeV. For example, for a typical gravitino mass ∼ 10 TeV, the BBN bound on the reheat temperature, T R (1 − 2) × 10 9 GeV, is consistent with the predictions displayed in figs. 2-4. Therefore, an unstable long-lived gravitino is viable for a wide range of reheat temperature described above.
For an unstable short lived gravitino, m 3/2 > 25 TeV, the BBN bounds on T R are no more applicable. However, there is another constraint coming from the decay of gravitinos to the lightest sypersymmetric particle (LSP) χ 0 1 . In this case, the LSP relic density is given by, [48] Ωχ0 where m χ 0 1 is the LSP mass, and Y 3/2 is the gravitino yield given as, Requiring the LSP neutralino density to not exceed the DM relic density, eq. (35) gives an upper bound on the neutralino mass, which is consistent with the lower limit on the neutralino mass mχ1 0 18 GeV [50]. Using eq. (37), the predicted range of T R with successful leptogenesis (i.e., 10 10 GeV T R 2 × 10 7 GeV) can now be translated into a viable mass range for the LSP neutralino, Thus, an unstable short-lived gravitino scenario is feasible for a wide range of LSP neutralino mass. This is in contrast to earlier studies of µ-hybrid inflation with M S ∼ 0 [19,20] where the feasibility of this scenario relies on a non-minimal Kähler potential.

GUT Embedding and Metastable Cosmic Strings
In a SUSY SO(10) framework, a metastable cosmic string network could arise via the symmetry breaking chain, where SM ≡ SU (3) c ×SU (2) L ×U (1) Y is the SM gauge group. The first breaking yields monopoles carrying SM and U (1) χ magnetic charges. The second breaking yields CSs with the string tension determined by the U (1) χ symmetry breaking scale [23,28]. Assuming that the monopoles are inflated away, the resulting string network is effectively metastable and yields a stochastic gravitational wave background spectrum that we explore in section 10. Note that the χ charge coincides with B − L for Y = 0. For a recent discussion of metastable CS netwrok formation in other gauge groups see [51]. The metastable string network decays via the Schwinger production of monopole-antimonopole pairs with a rate per string unit length of [23,[52][53][54], where m ∼ M G is the monopole mass and κ cs quantifies the metastability of CSs network with √ κ cs ∼ 10 being the stability limit as the lifetime of CSs becomes larger than the age of the Universe. This parameter plays an important role in making predictions for the current and future GW experimental tests.

Stochastic Gravitational Wave Background from Metastable Cosmic Strings
The stochastic gravitational wave background (SGWB) emitted from the CS network is calculated in terms of the fractional energy density in GWs per logarithmic interval of frequency [55], Here, H 0 = 100h km/s/Mpc is the Hubble parameter today with h = 0.68 [10], and P n is the power spectrum of GWs emitted by the n th harmonic of a CS loop. Our predictions are based on the Blanco-Pillado-Olum-Shlaer (BOS) model [55,56] and cusps as the main source of GWs with P n given by [57], where Γ 50 is a numerical factor specifying the CSs decay rate and ζ is the Riemann zeta function. It is convenient to work in terms of redshift z with 1 + z ≡ a 0 /a(t) written in terms of the scale factor a(t) and its present value, a 0 . The number of loops emitting GWs, observed at a given frequency f is defined as [55], The integration range corresponds to the life time of CSs network, from it's formation at z max 1 until it's 1 T R is taken to be around 10 9 GeV. (1+z)f . The loop density is defined by considering their formation and decay in different epochs. In a radiation dominated era the loop density is given by [57], with ≤ 0.1 t, whereas in the matter era it is given by [57], with < 0.18t. Lastly, the number density of loops produced during the radiation era, but radiating during the matter era is given by [57], with < 0.09t eq − Γ Gµ cs t.
The cosmological time as a function of z is written as [55], and the Hubble rate at redshift z is given by [55], where Ω m,0 = 0.31, Ω r,0 = 4.15×10 −5 h 2 and Ω Λ,0 = 1−Ω m,0 are the present values of matter, radiation and dark energy densities respectively, obtained from a standard flat ΛCDM model [10]. The function (z) defines the change in the expansion rate of the Universe due to annihilation of relativistic species at earlier times and is given as [59], where g * (z), g S (z) are the effective numbers of relativistic and entropic degrees of freedom respectively, at redshift z, and g * (0), g S (0) are their present values. The evolution of the effective degrees of freedom with redshift is shown in figs. 5a and 5b, both for the SM and MSSM 2 .
Here the parameter κ cs increases from left to right with √ κ cs = 10 representing the stable CS limit and Gµ cs decreases from top to bottom. The NANOGrav/PPTA [24,65] 2σ (1σ) posterior contours are shown by solid (dot-dashed) brown/orange region with a broken power law fit. The gray window in the upper right corner is the region excluded by the CMB constraint and is only applicable to CSs with lifetime longer than CMB decoupling, FIGURE 7: 7a shows the GW spectra explaining the NANOGrav excess at 2σ for allowed values of Gµ cs with √ κ cs = 7.9 to 9 and 10 as indicated with corresponding color ♠ markers in fig. 6 i.e., √ κ cs ≥ 8.6. The pink shaded region representing successful leptogenesis, DM and inflation is congruous with the gray bounded region of figs. 2-4. In fig. 7a the spectra of GW are shown for values of √ κ cs and Gµ cs which lie within the 2σ region of NANOGrav as indicated with corresponding color ♠ markers in fig. 6. Ignoring dependence on the effective degrees of freedom, the behavior Ω GW ∝ f 3/2 (Ω GW ∝ f 0 ) is achieved at the low (high) frequency range for the GW spectrum. This corresponds to the range γ cs ∼ 3.5 − 5 via eq. (51) while explaining most of the predicted region shown in fig. 6. The CS loops produced during the matter era and loops produced during the radiation era, but radiating during the matter era, become somewhat important in the low frequency range for γ cs 5, but this region lies outside of the NANOGrav bounds. Due to pair production of GUT monopoles, the metastable long strings on superhorizon scales and string loops and segments on subhorizon scales cause a SGWB which we have not considered here but can be seen in [66].
For a detailed comparison of the various present (solid) and future (dashed) experiments, the allowed values of √ κ cs and Gµ cs are shown in fig. 7b. It is interesting to note that metastable CS with √ κ cs ∼ 8 − 9 allow Gµ cs ∼ 10 −9 − 10 −6 . Therefore, the gravitino DM scenario with successful leptogenesis in µ-hybrid inflation, combined with NANOGrav SGWB, leads to the predicted range M G ∼ 10M ∼ 10 16 GeV − 10 17 GeV for Gµ cs ∼ 10 −9 − 10 −6 . However, at larger frequencies, the range Gµ cs ∼ 10 −8 − 10 −6 is in some tension with the latest bounds from LIGO O3 [41]. But this tension still involves some theoretical and experimental uncertainties. A non-standard thermal history [57,[67][68][69][70], or late production of the CSs [71] can ameliorate this tension. Besides NANOGrav, we also show in fig. 7b the observable region lying within the sensitivity bounds of future experiments, such as Einstein Telescope [63] and LISA [62].

Conclusion
We have explored the µ-hybrid inflation in a U (1) B−L extension of the MSSM by considering both the linear and quadratic soft SUSY breaking terms with special focus on the parameter space described by |M S | m 3/2 . A wide range of the gauge symmetry breaking scale, 6 × 10 14 M/GeV 10 16 is predicted for successful non-thermal leptogenesis and stable gravitino as a viable dark matter candidate. This parameter range corresponds to a stochastic gravitational wave background from a metastable cosmic string network with tension Gµ cs ∼ 10 −9 − 10 −6 . Such a metastable cosmic string network can arise in a grand unified theory with U (1) B−L embedded in an SO (10) model. An order of magnitude splitting is predicted between the GUT and B − L breaking scales. This connection certainly provides a unique opportunity to probe the seesaw mechanism and leptogenesis with gravitational waves [72][73][74][75]. The metastable cosmic string network lies within the 2σ NANOGrav 12.5 year data/PPTA and is also within reach of future gravitational wave experiments.