Dark Fermions and Spontaneous $CP$ violation in $SU(2)$-axion Inflation

Remarkably, if $CP$ was spontaneously broken in the physics of inflation, fermions would notice and remember it. Based on that, we present a new (non-thermal) mechanism for generating self-interacting dark Dirac fermions prior to the Hot Big Bang. The non-Abelian gauge fields and axions are well-motivated matter contents for the particle physics of inflation. In this background, we analytical study Dirac fermion doublets charged under the $SU(2)$ gauge field and use point-splitting technique to regularize the currents. We show that the non-trivial $CP$-violating vacuum structure of $SU(2)$-axion models naturally leads to an efficient mechanism for generating massive fermions during inflation. The size of the fermionic backreaction and the density fraction of dark fermions put upper bounds on the fermion's mass. For a GUT scale inflation, the generated dark fermions, only gravitationally coupled to the visible sector, can be as heavy as $m\lesssim 10 TeV$.


Contents
1 Introduction 2 2 Setup 7 2.1 Ψ + and Ψ − decomposition 9 3 Non-trivial vacuum structure and fermions 12 3.1 C-symmetry, P-and CP-violation 1 Introduction The current cosmological observations are in significant agreement with the general concept of inflation paradigm [1][2][3][4] and several of its key predictions have been already confirmed by the cosmic microwave background (CMB) and large scale structure (LSS) data [5]. Besides, the upcoming ambitious missions, e.g., LiteBIRD [6,7] and CMB Stage-4 [8], will soon provide us with even more data about the early Universe. However, our theoretical understanding of the particle content of inflation is still vague and incomplete in that sense. Given the fact that the energy scale of inflation can be as high as the GUT scale to a few TeV, it might provide an outstanding opportunity to explore the high energy physics beyond the standard model (BSM) and possibly breakthrough discoveries. In the past decade, several seminal and groundbreaking works have been done to systematically extract more information from the cosmological probes about the particle content of inflation through n-point functions, i.e. cosmological collider physics [9][10][11][12][13][14][15]. In these studies, the source of the particle production is either the expansion of the Universe or the decay of massive particles, while (up to this point) the vacuum is assumed to be topologically trivial and CP -preserving. From the high-energy physics point of view, gauge fields and axions are very natural matter contents for high energy scales, e.g., the scale of inflation. Axion fields are abundant in theories BSM and hence well-motivated candidates for the inflaton field. The first model of axion inflation, natural inflation, has been proposed almost 30 years ago in [16]. Enjoying an almost shift symmetry, the axion effective potential is protected from dangerous quantum corrections which guaranteed the flatness of the potential. Among the well-motivated axion inflation models, one can mention axion monodromy [17,18]. See [19] for a review on axions in inflation, and see [20] for a review on axions in cosmology. Besides, some theoretical issues of de sitter vacuum suggest not only the naturalness but the necessity of axions [21,22].
Another natural matter content for inflation is non-Abelian gauge fields. Thanks to the isomorphy of the su(2) and so(3) algebras, any non-Abelian gauge field can acquire a homogeneous and isotropic vacuum field configuration (VEV) in its SU (2) subsector. Moreover, their effective theory can break the conformal symmetry of the Yang-Mills theory, e.g., by coupling with the axion (to prevent its a −4 decay) and let the gauge field enjoy a slow-roll dynamics with an almost constant energy density during inflation. Besides, the non-Abelian gauge field perturbed around its isotropic VEV has an extra spin-2 degree of freedom which is linearly coupled to the primordial gravitational waves. Such SU (2)-axion vacuums spontaneously break P and CP during inflation [23] (See Fig. 1 for an illustration) and predict a CP -violating pre-Hot Big Bang physics. Since the seminal paper [24], it is well known that broken discrete symmetries cause cosmological domain wall problem which can dominate the Universe if the symmetry is broken only spontaneously and happens after inflation [25]. Given the fact that in SU (2)-axion models, the scale of symmetry breaking is high, and prior to the Hot Big Bang, they avoid the domain wall problem. The theoretically well-motivated and phenomenologically rich setting for SU (2) inflation models have been discovered in [26,27] which showed that non-Abelian gauge fields can contribute to the physics of inflation and the first realization of this class of models has been introduced, i.e., Gauge-flation. The next realization with the same phenomenology has been introduced as Chromo-natural inflation [28]. In addition to the original models, several more inflationary models with the SU (2) VEV have been proposed and studied which share the same key features. For a review on gauge fields in inflation see [29] and for a recent overview on the SU (2)-inflation models so far in the literature and their classification see Sec. 2 of [30] and the references therein. Figure 1. The SU (2)-axion vacuum structure. For each given parameter set and energy density, there are two SU (2)-axion field configuration vacuums, Eq.s (2.7)-(2.8), which are related by the parity, i.e. slow-roll trajectories A and A P . The two have identical background cosmologies; however, vacuum spontaneously breaks P at the level of perturbations of the fields with spin, e.g., fermions and spin-2 fields. For more details about the vacuum see Sec. 2 of [30].

Slow-roll
Despite their differences in details, all of the inflationary models with the SU (2) VEV share several key features as robust consequences of having non-Abelian gauge fields in the physics of inflation. The illustration in Fig. 2 summarizes some of them. The distinctive features of the non-trivial VEV are as follows. 1) The perturbed gauge field has a new chiral tensorial mode (spin-2 field) 1 which is linearly coupled to the gravitational waves, and generates chiral primordial gravitational waves [26,31,32], parity-odd correlations of CMB, i.e., non-zero T B and EB [33], and a possible large tensor bispectrum [34,35]. In [30], the backreaction of the spin-2 field on the background fields has been studied analytically for all of the models so far in the literature. 2 2) Thanks to the SU (2)-axion vacuum which spontanously broke CP in inflation and the new spin-2 field, this class of models has a non-trivial topology with RR = 0 which through the gravitational anomaly in SM, ∇ µ J µ l = N L −N R 384π 2 RR, provides a natural leptogenesis mechanism during inflation [23,37,38]. 3 Assuming type-I seesaw mechanism; the generated lepton number is proportional to the energy density of the gauge field during inflation, which provides the source of CP violation and the mass of the heaviest right-handed neutrino [37]. 3

) It provides a new mechanism for
1 For details about the nature of this tensorial mode in the SU (2) gauge field, perturbed around its isotropic and homogenous background, see App. B.1 of [30]. 2 At second order in perturbation, the spin-2 field couples to the curvature perturbation and contributes to the scalar power spectrum and bispectrum. For the chromo-natural model this effect has been studied in [36]. 3 The idea of inflationary leptogenesis through the gravitational anomaly in SM was first proposed within the context of modified gravity in the seminal work [39]. Observation Figure 2. The SU (2)-axion inflation setup is a complete beyond standard model (BSM) which provides new mechanisms to explain and relate seemingly disparate phenomena in inflation and pre-Hot Big Bang (e.g. massive dark fermions, new spin-2 field and chiral gravity waves) with late time cosmological observations (dark matter, matter asymmetry, and primordial gravitational waves). Namely, the non-trivial CP -violating vacuum leads to distinctive features for the spin-2 and fermion fields coupled to it. The solid (dashed) back arrows show the particle production by (backreaction of the particles to) the vacuum. The colorful arrows show emerged phenomena due to the corresponding inflationary particle production. The dotted lines indicate unexplored directions.
the non-thermal generation of heavy, self-interacting dark fermions during inflation. The extensive analytical study of this pre-Hot Big Bang fermiogenesis mechanism is the focus of the current work. 4) The produced heavy dark fermions might decay to the visible sector though some axion-mediated and/or direct feebly interactions which given the spontaneous CP -violation of the vacuum can provide a new mechanism for baryo-and leptogenesis. This promising possibility is beyond the scope of the current work, and we leave for future study. It worth mentioning that the SU (2) VEV can also produce scalar (Higgs) doublets coupled to it [40]. For a comparison between the scalar and fermion production in this setup see Sec. 7.1.1.
In de Sitter, gravitational production of massive fermions has been studied in [41], their Schwinger production due to a constant and homogenous U (1) field is studied in [42], and massive fermions coupled to slowly evolving axionic backgrounds was studied in [43,44] and the gravitational wave power spectrum induced by that fermions is studied in [45]. Including the SU (2) VEV in the physics of inflation, SM lepton and baryon production through gravitational anomaly has been studied in [23,37] (see also [46][47][48]), while massless fermions (in the absence of the axion) has been considered in [49]. Recently in [50] by the author, a doublet of massive Dirac fermions in SU (2)-axion setting charged under the 〈A μ 〉 〈∂ μ φ 〉 Figure 3. The processes underlying inflationary fermion production studied in this work. The gauge field vacuum naturally induces fermionic currents through iΨ / DΨ ⊃ g AΨ γ µ A µ Ψ (top). Moreover, fermions can possibly be coupled to the axion by effective interaction ∂ µ ϕΨγ µ γ 5 Ψ (bottom). SU (2) gauge field has been introduced, and the fermion backreaction to the background gauge field J µ a , has been studied numerically for a small part of the parameter space, i.e. m ∼ H. The current work presents an extensive study of the setup by computing all of the non-vanishing fermionic induced currents analytically and shows that this setup naturally leads to a new (non-thermal) mechanism for generating massive fermions during inflation. We regularize the UV divergent momentum integrals by the point-splitting regularization skim.
The paper is structured as follows. In Sec. 2 we review the setup of a dark Dirac fermion charged under the gauge field in a generic SU (2)-axion inflation. In Sec. 3, we discuss the symmetry structure of the SU (2)-axion vacuum and its consequences on the fermionic sector. In Sec. 4, we solve for the mode functions in (quasi) de Sitter. Next in Sec. 5, we take a closer look at the solutions and investigate their generic features. Sec. 6 presents the analytical form of the fermionic currents after point-splitting regularization. Sec. 7 discusses the main results and compares the inflationary dark fermion production in different limits. In Sec 8, we take a quick view on the (non-thermal) dark fermion production mechanism naturally arises in our non-trivial vacuum. We finally conclude in ??. Technical details of the computations as well as the underlying mathematical tools are provided in App.s A-D. Readers interested in the main results may proceed directly to Sec.s 7 and 8.
Notations and conventions: In this work, we adopt the notation introduced in [50] as much as possible. Although in some cases we use a more convenient notation to better capture the novel aspects of the setup. Here we deal with 8, 4 and 2-component spinors which are acted upon by 8 × 8, 4 × 4 and 2 × 2 matrices, respectively. We place a tilde (∼) on top of 8 dimensional spinor and matrix. The 4-dimensional spinor and matrix remains unchanged, while the 2-dimensional ones are written in boldface. Moreover, I n represents the n × n identity matrix and the gamma matrices are in the Dirac representation unless otherwise stated. We use the Einstein summation notation, i.e. repeated indices (one upper and one lower) are summed. Greek letters starting from the middle of the alphabet, i.e. µ, ν, . . . , are used for the space-time indices, whereas the starting ones, i.e. α, β, . . . , present the indices of the tangent space (non-coordinate) bases. In particular, the tetrads are defined as g µν ≡ η αβ e µ α e ν β . As a result, γ µ = e µ α γ α . We recall that the FRW tetrads can be defined as e µ α = a −1 (τ )δ µ α and γ α are the flat space Dirac matrices which obey the Clifford algebra where η αβ is the Minkowski metric with signature (− + ++). The indices of γ α and σ i are lowered as γ α = η αβ γ β and σ i = δ ij σ j . The gamma matrices in the Weyl representation are presented as γ α W . The 4-dimensional Hermitian adjoint row spinorΨ is given asΨ ≡ Ψ † γ 0 . Besides, γ 5 = iγ 0 γ 1 γ 2 γ 3 . The 4 × 4 helicity operator, h(k α ), is defined as (1.1) T a are the generators of the SU (2) group in the fundamental representation [T a , T b ] = i c ab T c , and we call the SU (2) index color. The generators of the SU (2) and the Pauli matrices are related as T a = 1 2 τ a . As a mathematical tool to take care of the entanglement of the color and Lorentz indices in our setup, we define the color-spin helicity operator (c-helicity) as in which k α is the momentum of the mode with k = | k|. Both σ i and τ a are the Pauli matrices. However, σ i represents the spin operator which carries a Lorentz index while τ a represents the generators of the gauge group and acts on color indices. Whenever convenient, we also use τ i = δ i a τ a . Throughout this work we will use three different representations of the doublet fermions, i.e. isospin, Weyl, and c-helicity frames which are related by two unitary operators,T 1 and T 2 . We present these mappings in the following table.

Representation
Isospin Here, Ψ 1 and Ψ 2 are the eigenstates of the isopin, while Ψ L and Ψ R are chiral states which are the eigenstates of γ 5 , i.e. γ 5 Ψ L,R = ∓Ψ L,R . Finally, Ψ ± are the eigenstates of the c-helicity operator defined in Eq. (1.2). More precisely, we have h(k α )Ψ ± = ±Ψ ± . Moreover, the subscripts s and superscripts p are labels of the c-helicity frame. Throughout this work, the superscript ± denotes the plus/minus subspace, while the subscript ± denotes the corresponding 2d helicity polarization states, i.e.
besides, the subscripts C and P denote the charge conjugated and parity transformed quantities respectively. In constructing physical observables associated with Grassmann variables and computing the expectation values, we use the antinormal ordering on the creation and annihilation operators of the spinors, . . ., defined as Finally, the point-splitting 4-vector which splits each point in the space-time into a forward and backward point is denoted by ε µ = (ε, ε i ). The forward/backward (physical) coordinates corresponding to X µ and the splitting 4-vector ε µ , are X µ f / X µ b defined as and we will eventually take ε µ to zero. Notice that the slow-roll parameter is shown by = −Ḣ H 2 .

Setup
In this work, we consider an SU (2) doublet fermion in an SU (2)-axion inflation setting. This setup has been recently introduced in [50]. We consider a doublet of Dirac fermions charged under the gauge field of the SU (2)-axion model with mass m and / D as For the sake of completeness, we also consider the possible effective interaction of the axion with the fermion as where ϕ is the axion field, f is the axion decay constant, λ is the dimensionless coefficient of the Chern-Simons interaction term of the axion, and β is a dimensionless parameter.
Here, ⊗ is the Kronecker product. For definition of ⊗ see App. A. Moreover, the quantity β can be of order unity or lower. Here, J µ 5 is the chiral current given as We assume slow-roll inflation with (quasi) de-Sitter metric where τ is the conformal time, H is the (almost) constant Hubble parameter during inflation, a(τ ) is the scale factor, and H ≡ aH.
Moreover, we assume a slowly-evolving homogeneous and isotropic SU (2) gauge field in the temporal gauge as [26] 4 and a homogeneous axion field where ψ(τ ) and ϕ(τ ) are two slow-varying pseudo-scalars during slow-roll inflation. Combination of the slow-varying SU (2) gauge field Eq. (2.7) and axion Eq. (2.8) make our SU (2)-axion vacuum (See Fig. 1). The class of inflationary models with such SU (2) VEV has several different realizations. See section 2 of [30] and the references therein for a recent comprehensive review on the models. In the perturbation sector, the SU (2)-axion vacuum appears in terms of the following dimensionless quantities 5 It is also useful to define a dimensionless parameter from the mass of the fermion as Moreover, it is more convenient to work with the canonically normalized field which is the comoving fermionic field. In the Fourier space, we can expandΨ(τ, x) as 12) and the fermion theory reads as (2.13) 4 The stability of this homogeneous-isotropic gauge field solution against initial stochastical anisotropies of the Bianchi cosmology has been studied in [51,52]. Interestingly, although such systems can in principle acquires anisotropic hair [53,54], the homogeneous-isotropic gauge field ansatz is the attractor solution (For the massive gauge field case see [55]). 5 Notice that ξϕ is different from ξ ≡ λ∂τ ϕ 2af H in the axion inflation backgrounds. In fact, ξ and ξϕ are related as ξϕ = βξ, where β can be of order unity.
The above theory is in terms of 8-spinors wherein the color and spin indices are entangled. The reason is the homogenous-isotropic gauge field's VEV in Eq. (2.7) in which the gauge group index associated with su(2) and spatial Lorentz index associated with so(3) has been identified to restore the rotational symmetry in the background (For a detailed discussion see [29]). The next step would be to find a frame in which the system is possibly reducible to lower dimension subspinors.

Ψ + and Ψ − decomposition
Up to this point, we were in the isospin frame. However, the ideal frame, would be the one in which the three of I 2 ⊗(γ 0 γ i ), I 2 ⊗k i γ i , and τ i ⊗γ i are diagonal. If exists, such a frame is made of two copies of the common eigenstates of I 2 ⊗k i .σ i , and τ i ⊗σ i . However, these two operators have only two common eigenstates. Thus, the theory is only block-diagonalizable and can be decomposed into two subspaces.
As a mathematical tool for any given 4-momentum k α , we define the color-spin helicity operator (c-helicity) as in which the first Pauli matrix is the su(2) generator while the second one is the spin operator. As we will show shortly, it is easier to expand the system in terms of the eigenstates of c-helicity, which automatically takes care of the entangled color and spin indices in the action (2.13). The orthonormal eigenstates of c-helicity satisfy in where p = ±1, s = ±1, and E p s (k α ) are (2.16) in which σ α andσ α are and their indices are lowered with the Minkowski metric, i.e. σ α = η αβ σ β . In a similar way, we define τ α andτ α from τ i . Besides,ǩ α is a four vector given aš where k = √ k.k. Notice thatǩ α is the four momentum of the massless field, but for the massive cases it is just a mathematical tool. It is straightforward to see that E p s (k α ) are eigenstates of helicity operator as well and satisfy the orthonormality condition (2.20) The p = +1 elements are also eigenvectors of τ i ⊗ σ i in Eq. (B.6) However, the p = −1 elements does not satisfy eigenstate equations with τ i ⊗ σ i operator. Thus, the system is not fully, but only block diagonalizable. The E p s (k α ) 4-vectors make an orthonormal basis and two copies of them define a 8 × 8 frame, i.e c-helicity frame. As it is implies by the form of E p s (k α )s, the color and spin are totally mixed in this frame. It is shown in [50] by the author that the system is decomposed into two independent subsystems. For self-sufficiency, we provide the details in App. B.1 and here we only present the final forms. The plus subspace is made of two copies of E + s (k α ), and the minus subspace is made of two copies of E − s (k α ), i.e.
such that the theory is given as The system is reducible into two 4-Dirac spinors, Ψ + and Ψ − , thanks to the isospin symmetry as well as the isotropic-homogeneous configuration of the gauge field's VEV which identifies the su(2) gauge index and so(3) rotation index. In particular, in case that the mass term breaks the isospin symmetry, i.e., the isospinors have different masses, the Ψ + k and Ψ − k fields would be coupled by terms proportional to the difference of the isospinor masses. Hence the reduction is not possible in the absence of the isospin symmetry.
Our Dirac fields can be expanded as where ψ ±↑ s (τ, k) and ψ ±↓ s (τ, k) are mode functions and E s with s = ±1 are the two-spinor polarization states Note 1: the superscript ± denotes the plus/minus subspace, while the subscript ± denotes the corresponding 2d helicity polarization states, E ± . Note 2: being already in the helicity frame for the given momentum, the 2-spinor polarization states are k-independent. As a result, the plus spinor is decoupled in terms of the polarization spinor E s and it is makes two pairs of coupled field equations for each polarization. On the other hand, the minus spinor is not diagonalizable due to the extra (time dependent) term proportional to γ 1 in the action (2.24). As a result, the minus subsystem includes four coupled field equations. In the limit of either well inside the horizon, i.e. k H, or ξ A /µ m 1, this term is negligible and the minus system approximately decoupled into two pairs of coupled field equations.
Our fermionic sector generates the fermion current and isospin current as as well as axial current, J µ 5 , in Eq 2.28 and axial isospin current Among them, J µa is the Noether current which backreacts on the background equation of the gauge field as 7 where the left hand side is the background field equation of the SU (2) gauge field and L inf is the specific theory of the SU (2)-axion inflation model. Moreover, the axial current backreacts on the axion background equation as where the left hand side is the background field equation of the axion, and the ∇ µ J µ 5 is given as In the right-hand side, the first term is the tree-level effect of explicite breaking of the chiral symmetry by the mass term, while the second term is the loop effect due to the well-known chiral (Adler-Bell-Jackiw) anomaly [56,57].

Non-trivial vacuum structure and fermions
To gain more insight on the nature of our dark fermions and capture some of the novel aspects of the vacuum, here we explore the action of the discrete symmetries C, P and CP , as well as the continuous chiral symmetry on fermions. Here, just by using the symmetry structure of the system, we extract the general features of the observable quantities. Later in Sec. 6, we explicitly compute these quantities and confirm this qualitative analysis. The chiral charge and the fermion backreaction to the background gauge field equation, Eq.s (2.28) and (2.27), in the Weyl frame can be written as 8 in which the R/L subscription denotes the contribution of the right-/left-handed fields to that quantity. The net fermion number Eq (2.29) is given as and the axial isospin current can be written as

C-symmetry, P-and CP-violation
As showed in Eq. C.1, our fermion theory in Eq. (2.13) is invariant under the action of the charge conjugation operator,C, defined as In other words, bothΨ, and its charge conjugated field Ψ c ≡CΨ, obey the same theory (See Eq. (C.8)) and the theory is C-symmetric. The 4 × 4 operator C is the charge conjugation operator of each plus and minus subsystems. However, parity is different since both ϕ (axion) and ψ (the effective field value of the SU (2) gauge field) in Eq.s (2.7) and (2.8), are pseudo-scalars. As a result, parity is spontaneously broken in the action (2.13). Under the action of parity, this theory goes to its twin one with (see Fig. 1) As a result, both P and CP are spontaneously broken by the vacuum of the SU (2)-axion, i.e., the vacuum is CP -violating, and it has a non-trivial topology with non-vanishingRR 8 Notice that the chiral charge is defined in terms of J t 5 = ∂t ∂τ J 0 5 which justifies the extra factor of a inside the integral Eq. (3.1). and FF [23,37]. Violation of CP is the essential condition in any matter-antimatter asymmetric scenario. Then the out of thermal equilibrium state which is guaranteed by inflation ensures that the production of matter by one mechanism is not immediately compensated by its disappearance through the inverse reaction. This robust aspect of SU (2)-axion models make them natural settings for inflationary leptogenesis [23,37]. In the current fermionic sector with C-symmetry and CP -violation, we expect a vanishing net fermion number, J 0 = 0 and a non-zero chiral charge, J 0 5 = 0 which will be confirmed by the exact computation in section 6. There we will show that both J and Q 5 are directly proportional to the sources of the CP -violation, i.e., ξ A and ξ ϕ .

Chiral anomaly and conserved currents
In the massless limit, the fermionic sector enjoys the chiral symmetry in the classical level under the unitary chiral transformatioñ where t a is the isospin matrix (equals to T a for isospin doublet and I 2 for isospin singlet), and α a (x) is the parameter of the transformation. Namely, all four fermionic currents are classically conserved under the above unitary transformation. The fermion vector current and isospin current remain conserved at the quantum level as well In particular, J µa is the Noether current, and its conservation is the fundamental consistency condition for the gauge theory [58]. However, the axial vector and axial isospin currents receives corrections according to the Adler-Bell-Jackiw anomaly [56,57] respectively as and Therefore, in this setup the axial isospin current is not affected by the chiral anomaly and remains conserved. 9 In the zero mass limit, the only fermionic effect is the Abelian chiral anomaly in Eq. (3.9),

Fermions in SU(2)-axion inflation
Now we turn to solve the field equations and find Ψ ± k specified by the actions in the Eq.s (2.23) and (2.24). The Ψ p k spinors with p = ±1 can be decomposed as where a p s,k and b p s,k are creation and annihilation operators obey Note that the superscript p = ± denotes the plus/minus subspace, while the subscript s = ± denotes the corresponding 2d helicity polarization states, E ± . Moreover, the charge conjugation relates the field to its charge conjugated field as (see Sec. C.1) The 4-spinor fields, Ψ p k , are govern by first order differential equations. Therefore, they are each made of four linearly independent solutions. We will solve them by setting the Bunch-Davies vacuum as the initial condition for the fields Here, we define the rescaled physical momentum as Since the plus and minus spinors are decoupled, we study them individually.

Ψ + fermions
As we discussed in Sec. 2.1, the Ψ + k spinors are defined to be diagonalizable in terms of the c-helicity. The Ψ + k spinors are described by the linear differential equation specified by the action (2.23). Therefore, U + s,k (τ ) and V + s,k (τ ) in Eq. (4.1) can be decomposed in terms of the 2d helicity polarization states, E ± , defined in Eq. (2.26) as .
Moreover, using Eq. (4.4), we can read V + s,−k (τ ) as . (4.8) Using (4.7) in the action (2.23), we arrive at two sets of coupled field equations for u ↑ s and u ↓ s as The above coupled set of equations can be decoupled in terms of the following decomposition The coupled set of first order differential equations (4.9) and (4.10) can be decoupled into two second order differential equations for Y s and Z s as and µ + is The general solutions for Eq.s (4.12) and (4.13) are W and M Whittaker functions. Setting the Bunch-Davies vacuum as the initial condition for u ↑ s and u ↓ s and using Eq.s (A.9) and (A.10), we find Y s and Z s as (4.16)

Ψ − fermions
Next, we turn to Ψ − k modes, described by the action (2.24). The theory of the minus subspace is not diagonalizable in a time independent frame. Hence, the helicity eigenstates are only the eigenstates of the Lagrangian in the asymptotic past limit. Therefore, we use the ansatz where E s are the 2-spinor polarization states defined in Eq. (2.26). Since in the asymptotic past limit the system is diagonalized in this particular basis, we have Moreover, using the charge conjugation relation, we can read V − s,−k (τ ) as .

(4.19)
Using the ansatz Eq. (4.17) in the action (2.24), we arrive at the field equations below Note that for each given s and s , the equations with s = + and s = − are coupled. We make the decompositions (4.23) We can reduce the above coupled first order equations to pairs of coupled second order equations Unlike the plus spinors, the set of equations in minus system is still coupled and cannot be solved analytically. However, since the coupling term becomes relevant around and after horizon crossing, it is possible to solve the system analytically in the limit that the coupling is negligible, i.e. kτ 1 and/or ξ A /µ m 1. More precisely, we can decompose the fields as 10 Y s,+ (τ ) = y s, where Y s , and Z s are analytically solvable while the other four functions, y s,p (τ ) and z s,p (τ ), should be solved numerically. In the asymptotic past, we have and z s,− (τ ) up to a phase, .33)). The fields Y s and Z s satisfy (4.32) The above differential equations can be solved analytically in terms of Whittaker functions and setting the Bunch-Davies vacuum as the initial condition for u ↑ s,+ and u ↓ s,+ , Y s and Z s are given as Next, using the ansatz Eq. (4.27), as well as the first order differential equations of Y s and Z s in Eq. (4.23), we find a set of coupled first order differential equations for y s,s and z s,s . Given the analytical forms of Y s and Z s in Eq. (4.33) and the initial condition Eq.s (4.28)-(4.31), we solve the above coupled linear differential equations numerically. The results are presented in Sec. 5.2 ( Fig.s 8-15). In order to have a better insightful vision of the minus subsystem, next we discuss the qualitative behavior of Ψ − k and the y s,s and z s,s functions.

Geometry of the minus subspace
Unlike the plus subspace, helicity is not a good quantum number in the minus subspace. That is because of the off-diagonal term proportional to ξ A in its action in Eq. (2.24). Therefore, although both Ψ + and Ψ − spinors have similar asymptotic past limits, but the more minus spinors are stretched and approach the horizon, the more they experience the off diagonal time-dependent term. As a result, a minus fermion with a given initial polarization state E s , gradually rotates in the helicity space and turns to a mixed state of E s and E −s , i.e.
where the time dependent polarization states, E Y s (τ ) and E Z s (τ ), can be read off in terms of the helicity eigenstates, E ± , as All E Y s (τ ) and E Z s (τ ) polarization states start from E s at asymptotic past. Then, around and at horizon crossing, they make most of their rotation in the E + − E − plane. Here, for more insight, Fig. 4 shows schematic rotation of E Y Before going any further and numerically solve the minus subsystem in Sec. 5.2 ( Fig.s  8-15), let us discuss the qualitative behavior of these functions.
• The y s,− (τ ) starts from zero at asymptotic past and after a quick transition phase increases to a constant number less than one after horizon crossing, i.e. y 0s− . Its behavior can be roughly described by 1 2 y 0s− − tanh[α sτ −τ 0s ] + 1 with y 0s− 1.
• The z s,− (τ ) starts from one at asymptotic past and after horizon crossing approach a constant value either less or more (but of the order of) one. Yet it also follows a tanh-type behavior.

Mode functions, a closer look
In this part, first, we discuss the generic features of the analytical solutions Y s and Z s , which are in common between the plus and minus subsystems. Next, we present the numerical study of the four functions, y s,p and z s,p , which we introduced in the minus subspace to take care of the rotation of its mode functions in the helicity plane.

Ψ ± spinors: Y s and Z s
Here, we consider Y s and Z s mode functions which are in common between the plus and minus sub-sectors. In fact, in both subsystems, µ is pure imaginary, i.e. (see Eq.s (4.15) and (4.26)) while κ ± s andκ ± s in Eq.s (4.14) and (4.25), can be written as where κ ± I is the absolute value of the imaginary parts of κ ± s (andκ ± s ) and µ ± can be decomposed as Here, |µ ± |, m ± eff , and κ ± I are real parameters in terms of ξ A , ξ ϕ , and m H . As a result, we can present the analytical solutions in both subspaces, Eq.s (5.5) and (4.33), in a unified form and in terms of three independent parameters, ξ A , κ I , and m eff as Zs √ x with s = ± in terms ofτ for several values of parameters. Moreover, Fig. 7 shows the same quantities as well as in the large mass limit. Here we summarize the generic behavior of these mode functions.
• The frequency of these oscillations as well as the amplitude of Y ± increase with the increase of κ I . However, the amplitude of Z ± decreases with κ I .
• Increasing m eff without changing µ m , which is possible in the minus subspace, decreases Y − and Z − , while has negligible effect on Y + and Z + .
• Increasing µ m has negligible effect on Y s while it linearly increases Z s , i.e. Z s ∝ µ m .
• In the large mass limit, µ m ξ φ , ξ A and for each s = ±, Z s increases to a constant value very close to (but less than) Y −s . See the right panels of Fig. 7.
• For the Ψ − k , one needs y s,s (τ ) and z s,s (τ ) functions in addition to the Y s and Z s to fully identify this subsystem. That is the task of Sec. 5.2.

Ψ − spinors: y s,s and z s,s
In Eq (4.27) we decomposed the fermion mode functions of the minus subspace in terms of the analytical functions, Y s and Z s , as well as four unknown functions, y s,± and z s,± which should be solved numerically. In this part, we present the numerical solutions of these functions in terms ofτ and for several values of µ m , ξ A , and ξ ϕ parameters. Moreover, for later use, here we also present |y s,+ | 2 + |y −s,− | 2 , and |z s, Later we will see that these are the combinations which appearers in the currents of the minus spinors.
• Neutral fermions: In this limit the fermion is not charged under the SU (2) gauge field, and is only derivatively coupled to the axion. In Eq.s (4.27) which formulate the deviation of the minus spinors from the plus one, we have y s,+ (τ ) = z s,+ = 1, and y s,− (τ ) = ξ A µm z s,− (τ ) = 0. For the fermions which are coupled to the gauge field through ξ A , this limit is the same as going to ξ A = 0. See the ξ A = 0 (gray) curves in Fig.s 8-15. Therefore, the minus subspace is simply another copy of the plus subspace.
• Large mass limit: Two different regimes with large mass limit are presented in The generic features of the system at this limit is as follows. The y s,+ (τ ), and z s,+ (τ ) (y s,− (τ )) start from one (zero) at asymptotic past and after a fast transition phase each reduces (increases) to another constant number less than one after horizon crossing. The z s,− (τ ) starts from one at asymptotic past and after horizon crossing approach a constant value of the order of µm ξ A . All of the curves roughly follows a tanh-type behavior. Furthermore, in the very massive fermion limit, the combinations |y s,+ | 2 + |y −s,− | 2 , and |Z s, |Z −s,− | 2 (s = ±) start from one at asymptotic past and remain very close to it throughout the mode's evolution. • SU(2) Schwinger process: The system with no axion coupling, ξ ϕ = 0, is presented in Fig.s 14 and 15. This limit is simply the SU (2) Schwinger process during a (quasi) de-Sitter expansion.

Fermionic currents: point-splitting regularization
Working out the fermion fields in the previous section, we are now ready to compute the induced currents. The three currents induced by our fermions are presented in Eq.s (2.28) and (2.27). However, those expressions for currents are in their classical forms. Given the Graussman nature of the fermion fields, i.e.Ψ † α andΨ β do not commute, the quantum currents are given by the anti-symmetrization of the fermionic fields. Therefore, the physical fermion (vector) current is given as [59] J µ = 1 a δ µ α 0 [Ψ, I 2 ⊗ γ αΨ ] 0 .    It is more convenient to replace the commutation on operators with the antinormal ordering on the creation and annihilation operators, . . ., defined as Therefore, the physical fermion current can be written as Similarly, we can find the physical isospin current and chiral current from their classical forms in Eq.s (2.27) and (2.28).

Fermion number
Here we compute the net fermion number, J 0 . From Eq. (6.3), we find the contribution of the plus fermions to the net fermion number as where in the last equality we used Eq. (4.8). Similarly, the contribution of the minus fermions to the net fermion number is As a result, the net fermion number vanishes which is a consequence of the C-symmetry in the fermionic sector.

Isospin current and chiral charge
At this part, we compute the isospin current, J µa , and the chiral current, J µ 5 . The corresponding momentum integrals are UV divergent and need to be renormalized for which we will use the point-splitting regularization skim. The gauge invariant isospin and chiral currents with symmetric point separation are respectively given as and in which W(y, z) is the Wilson line defined as which is required to maintain the gauge invariance and trace is over the fundamental representation of the gauge group, P stands for path ordering, and x µ f and x µ b are the forward and backward point-splitting points, corresponding to the point-splitting physical points 11 respectively and ε µ is a constant 4-vector which we will eventually take to zero. Using the slow-roll inflation relations, Ha(τ ) − 1 τ and a(τ ) e Ht , we can read off the forward and backward point-splitting conformal times as and the spatial comoving coordinates as In terms of the comoving fields and using the plus and minus decomposition, the chiral charge in Eq. (3.1) can be written as where Q + 5 and Q − 5 are the contribution of the plus and minus fermions to Q 5 respectively. Similarly, the fermionic backreaction to the gauge field's equation in Eq. (2.29) can be written as where J ± are the contributions of the plus and minus spinors to J . Taking the ε µ goes to zero limit and using Eq. (B.31), Q ± 5 (x ν ; ε) can be written as . . . 0 . (6.18) 11 It is noteworthy to mention that in the expanding universe, the point-splitting, i.e.
should be done to the physical coordinates and not the comoving ones. Otherwise, the UV divergent terms would be powers of τ ε , not 1 ε . Such a mistake leads to an unphysical term of the form ln(a) in the currents which fakes an IR divergence. 12 More precisely, the explicate form of J − (x ν ; ε) is However, from the combination of the charge conjugation Eq. (4.4), as well as antinormal ordering in Eq.
It is noteworthy to mention that the first term in Eq. (6.15) is the contribution of ε i terms in the presence of the gauge field and the Wilson line. This term is the zeroth element of the Chern-Simons current and the result of the Abelian Adler-Bell-Jackiw anomaly [56,57] (See Eq. 3.9). In the momentum integrals on the right hand sides of Eq. (6.15) and (6.18), all the divergent terms are given in terms of the temporal separation. That is due to the spatial isotropy of our setup. Therefore, from now on, we set ε i to zero and only keep the separation in time parameter, i.e. ε. Interestingly, both the isospin backreactions, J ± , and chiral charges, Q ± 5 , are specified in term of the following dimensionless quantity More precisely, the currents can be written as Besides, the fermionic backreaction to the axion background in Eq. (2.30) can be written as In the following, we compute K + (x µ ; ε) and K − (x µ ; ε) which are the contributions of the plus and minus fermions respectively.

Ψ + fermions:
In the plus subspace and after using Eq.s (4.8) and (6.19), K + can be written as

(6.22)
This momentum integral is UV divergent and we use the point-splitting skim to regularize it. Here, we present the final results and the details of the computations are presented in App. D. We find that K + (τ, ε) has a log divergent term and can be written as where K + reg (τ ; ε) is the regularized current, given as in which ψ (0) (z) ≡ dΓ(z) dz is the digamma function and κ + I is the imaginary part of κ + , i.e. κ + I = 2ξ ϕ − 1 2 ξ A , and |µ + | = [µ 2 m + (κ + I ) 2 ] 1 2 . Notice that K + reg is proportional to κ + I and the quantity inside the curly brackets is even under parity.

Ψ − fermions:
For the minus subspace and after using Eq. (4.19), we obtain K − (τ, ε) as Since we only have y s,s (τ ) and z s,s (τ ) functions numerically, we can not work out the current analytically. However, as it is expected from Eq.s (4.23) and is showed in Sec. 5.2, in the limit ξ A /µ m 1, the combinations |y s,+ (τ )| 2 + |y −s,− (τ )| 2 , and |z s, |z −s,− (τ )| 2 are always close to one. As a result, we can approximate K − (τ ; ε) in this limit as Following the corresponding analysis in the plus subspace, we find the regularized current in the minus subspace, K − reg , as 2 . As we mentioned earlier, is only valid in the limit ξ A /µ m 1.
To summarize, we worked out the induced currents analytically by using the pointsplitting regularization skim. In the minus subspace, our analytical solution for J − and Q − 5 (given in terms of K − ) are approximations which are valid in ξ A /µ m 1 limit. To find K − reg (τ ) in other parts of the parameter space, one needs a full numerical analysis, e.g. [50]. The two approaches are only in agreement when the source of the particle production, κ ± I is large. The reason is in the use of a hard cut-off in the numerical computation of the momentum integrals carried out in [50]. As a result, the numerical approach subtracts the divergences correctly but misses the finite terms which appear by the regularization process.

Dirac fermions in (quasi) de Sitter
In the previous section, we worked out the analytic form of the fermionic currents sourced by an SU (2) gauge field and possibly derivatively coupled to the axion field in de Sitter. In this section, we take a closer look at the results and make a detailed comparison with other cases of fermions generated in different setups in de Sitter, i.e., neutral fermions coupled to axion, and Abelian fermions. We also compare the non-Abelian fermion Schwinger effect with its scalar counterpart.

Charged SU (2) fermions
This setup of Dirac fermion doublet is reducible into two subspaces of Dirac fermions in the c-helicty frame (Sec. 2)Ψ where the supercript +/− denotes the plus/minus subspace, and Ψ ± k are two 4-spinors. In Sec. 6, we showed that our non-Abelian fermions have vanishing fermion number, i.e.
as a consequence of C-symmetry in our setup. Then, we showed that the isospin current backreacts to the (background) SU (2) gauge field, J ≡ a 3 δ i a J ai = J + + J − , and also has a non-zero chiral charge, Q 5 ≡ aJ 0 5 = Q + 5 + Q − 5 . We regularized both of these quantities in Sec. D and found where K ± reg is a momentum integral which we computed analytically for the plus subspace in Eq.s (6.24) (and for minus subspace in ξ A /µ m 1 in Eq. (6.28)). During slow-roll, they are almost constant quantities in terms of the (slow-varying and CP -breaking) vacuum quantities ξ A , ξ ϕ in Eq.s (2.9), as well as the mass of the fermion µ m = m H . Thus, the backreaction to the axion background, B ≡ βλ 2f ∇ µ J µ 5 , is given as The K ± reg are directly proportional to (and an odd function of) κ ± I = (2ξ ϕ ∓ ξ A /2) which is the source of CP -violation and increase with the mass of the fermion, µ m , which is the source of chiral symmetry breaking. Thus, K ± reg vanishes in the massless limit where the system has classical chiral symmetry. Moreover, in the absence of the SU (2) gauge field and the axion interactions which are the sources of CP -breaking, the fermion fields are unpolarized and all the currents vanish. That is in agreement with our qualitative discussion based on symmetries in Sec. 3. In the large and small mass limits K + reg is In the limit that κ + I µ m , the K + reg (τ ) takes the following asymptotic form In Fig. 16, we plot K ± reg for fermion fields coupled both to the gauge field and axion which ξ ϕ = 2. Figure 17 shows K ± reg for fermions coupled to the isotropic SU (2) gauge field but in the absence of the derivative coupling with the axion, i.e. β = 0 and ξ ϕ = 0.

Light fermions
In the massless limit, K ± reg = 0, and the only non-vanishing quantities are the the chiral charge 5) and the fermion backreaction to the background axion field Namely, the only non-zero effect is the Abelian chiral anomaly (See Sec. 3.2).

Heavy fermions
The mass of the fermion and the energy scale of inflation cannot exceed the decay constant of the axion, f , which is the U V cutoff of the theory. Apart from that, the size of the fermion backreactions on the background gauge field and axion put upper bounds on it. Note that unlike the non-perturbative tunneling effects in Minkowski space, our current analysis is perturbative and relays on assuming a classical background with infinite energy that sources the fermions. That approximation is broken once the fermion backreaction becomes large. Following [30], here we constrain the parameter space of the system based on the size of the fermion backreaction. I) Backreaction to the background SU (2) gauge field: The validity of the perturbation theory and slow-roll dynamics requires that the backreaction to the gauge field should utmost be as small as slow-roll suppressed terms. Demanding that condition, we find which in the large mass limit, and requiring the RHS be smaller than 10 −2 H 2 ψ, gives H, (7.8) which for typical values of κ + I and ξ A of order one and g A ∼ 10 −2 − 10 −3 gives m 100H. II) Backreaction to the background axion field: Demanding that the femrion backreaction to the axion field equation should utmost be as small as slow-roll suppressed terms implies In SU (2)-axion models during slow-roll inflation, we have ξ A ∼ ξ [30] which using in the above, it gives For instance, for f /H ∼ 10 2 we have ξ A ∼ m H < 600.

Charged scalars vs charged fermions
Here we compare the scalar and fermion backreactions to the background SU (2) gauge field. As for the backreaction of the charged scalars, we use the analytical formula worked out in [40]. In Fig. 18, we show the regularized scalar and fermion backreactions. Interestingly, our vacuum produces more fermions than scalars in i) large ξ A and ii) large mass µ m > 1, i.e. m > H limits. More precisely, the asymptotic form of J reg in the large ξ A for the scalars decreases like 1/ξ A while for the Dirac fermions increase as ξ 2 A . That is because the SU (2) gauge field VEV adds an additional mass term to the spin-0 field in such a way that the charged scalar has negligible particle production in the strong field limit [40]. Moreover, by the increase of their mass, the scalar current decreases like 1/µ m while the fermionic one increases like µ 2 m . As a result, the SU (2)-axion inflationary models can efficiently produce massive fermions.   Figure 18. Comparing the scalar and fermion backreactions to the background SU (2) gauge field. The solid lines denotes J + reg for charged Dirac fermions (with ξ ϕ = 0) and the dotted lines represent the corresponding backreaction of the charged scalar fields worked out analytically in [40].

Neutral fermions coupled to axion
As a limit in our setup, here we focus on neutral Dirac fermions coupled derivative to the axion. 13 In the absence of the interaction with the gauge field, the minus subsystem is simply a second copy of the plus subsystem, hence K − reg = K + reg = K reg . Thus from Eq. (7.1), its backreaction to the gauge field vanishes, i.e.
Moreover, the fermion number is zero, i.e.
However, it still has a non-zero chiral charge while the backreaction to the axion field equation is given as where K reg is in which κ I = 2ξ ϕ and |µ| = (µ 2 m + κ 2 I ) 1 2 . The chiral charge increases with the mass of the fermion like µ 2 m and the axion parameter like ξ 2 ϕ . The K reg in this case is presented in Fig.  19.

Abelain fermions, small mass limit
For the sake of completeness, let us now consider the case of a Dirac fermion coupled to an Abelian gauge field in the small mass limit in de Sitter. The aim is to explore the relation of the mass of the fermion and its backreaction in the U (1) gauge field case qualitatively to explore if it also shows the increase with m behavior which we found for the other cases. Consider a U (1) gauge field in the temporal gauge, with a constant electric field as 13 Neutral Dirac fermions derivatively coupled to the axion, its backreaction to the axion ∇µJ µ 5 = −2imΨγ5Ψ, and the fermionic effects on the curvature perturbations have been extensively studied in [44]. The U V divergent momentum integrals are regularized using adiabatic regularization for fermions. The authors reported a singularity in the massless limit and used a change of basis to remove it. In this work, we did not find any singularities. In the massless limit in particular, we realize that the neutral fermions acquire conformal and chiral symmetry and all the observable currents vanish. and a Dirac fermion which is charged under it with the theory where Ψ is the canonical Dirac field (for which the spin connection disappears), and m is the mass of the fermion. In this case, the Noether current is 15) and the fermion backreaction to the U (1) gauge field is given as Moreover, in this case FF = 0, and hence both the vector and axial vector currents are conserved. Namely, the massless fermions have chiral symmetry, hence their backreaction to the gauge field should vanish. It implies that the backreaction should be directly related to the mass of the fermion which is the source of the chiral symmetry breaking 14 As a result, our simple argument implies that even in the U (1) gauge field case, the fermionic backreaction should increase with the mass of the fermion.

Fermion dark matter, a quick view
Up to this point, we analytically studied the (dark) Dirac fermions charged under the SU (2) VEV in a generic SU (2)-axion inflationary model. We realized that the spontaneous CPviolating vacuum structure of this class of models leads to a new and efficient mechanism to generate non-thermal massive dark fermions with a preferred helicity state, h-asymmetric, during inflation. Moreover, the generated dark fermions can be very massive, which depending on the efficiency of their possible feebly interactions with the visible and lighter dark sectors, may decay to leptons and baryons or other dark particles after inflation. In other words, these Dirac fermions may replace the right-handed neutrinos in the standard thermal leptogenesis scenario and hence provide a novel natural scenario to explain the observed matter asymmetry in the Universe. However, this possibility, as well as the late time macroscopic properties of the produced dark fermions as a dark matter sector, are all highly model-dependent. They depend on the value of the gauge coupling, g A , the possible mass of the gauge field at some symmetry breaking scale after inflation, and the details of the possible feebly interactions of the dark fermions with the visible and dark sectors which is beyond the scope of the current work. We leave these important questions for future exhaustive study. However, let us here roughly estimate the gist of the simplest possible scenario for fermio-genesis in this setup. The three necessary Sakharov conditions to generate an imbalance of baryonic matter are i) C-symmetry and CP-symmetry violation, ii) interactions out of thermal equilibrium, and iii) baryon number, B, violation [61]. One can generalize it for any form of matter by replacing the third condition by the violation of its corresponding fermion number. The fermion doublets charged under our SU (2) gauge field which we studied so far are C-symmetric. That, however, can be easily and naturally violated by the mass term, nonperturbative effects, and especially interactions with the SM fermions, which is beyond the scope of the current work. In the absence of these effects, the generated fermions are not asymmetric in the standard sense but fermions with h-asymmetry.
For simplicity, we assume that the dark sector is made of one generation of massive dark fermion doublets,Ψ. We also assume that they are only gravitationally coupled to the visible sector. Thus, dark fermions are stable. The dark fermion has a constant energy density during inflation where H is the Hubble parameter during inflation. Notice that the contribution of the ± spinors into ρ df (which is positive definite) should be directly proportional and an even function of κ ± I which is the source of the particle production. After the end of inflation, the source of fermion production (gauge field and axion) disappears and therefore ρ df decays like ordinary matter, i.e.
Here, in order to estimate the energy density of the reheating, we consider the phenomenological reheating model below where σ is the efficiency of the reheating process and relates ρ reh and the energy density at the end of inflation, ρ inf . From that, one can read the number density of photons during reheating as where g eff = 427/4 is the number of relativistic degrees of freedom at the time of reheating, and ζ(z) is the Riemann zeta function, ζ(3) = 1.2. Now, we define the parameter θ df as the ratio of the dark fermion energy density to the number density of photons as For a scale of inflation as high as GUT scales, H ∼ 10 −6 M Pl , and σ of order one, that gives m ∼ 10 T eV . Notice that it is much lower than the upper bound imposed by the size of the backreaction in Eq. (7.8) and higher than the electroweak scale, E EW ∼ 246 GeV .

Discussion
Models of inflation with an axion and SU (2) gauge field have very rich phenomenology. In this class of models, the SU (2) gauge field has a (slow-varying) homogeneous and isotropic VEV with an almost constant energy density during inflation. Both of the axion and effective VEV of the gauge field are pseudo-scalars (see Eq.s (2.7)-(2.8) and Fig. 1). As a result, both P and CP are spontaneously broken by the vacuum. Although the background evolution is insensitive to that symmetry breaking, cosmic perturbations with spin, e.g., fermions, and spin-2 fields notice that and become chiral (see Fig. 2). In particular, the Pviolation leads to (a short phase of tachyonic) production of the spin-2 field of the perturbed gauge field which produces chiral primordial gravitational waves. Fermions are even more sensitive to the structure of the vacuum as a consequence of the celebrated AtiyahâĂŞSinger index theorem. As first pointed out in [23,37,38], the spontaneous CP -violation and nontrivial topology of this vacuum with RR = 0, provides a natural leptogenesis setting via the gravitational anomaly in SM. In this work, we unveiled yet another unique aspect of such SU (2)-axion vacuum for fermions which might open a new window toward particle cosmology.
In the context of SU (2)-axion inflation models, we considered a (dark) massive Dirac fermion doublet charged under the SU (2) gauge field and (for the sake of completeness, possibly) derivatively coupled to the axion. This setup has been recently introduced in [50] in which the fermion backreaction to the gauge field background has been numerically computed for a small part of the parameter space, i.e., m ∼ H. Here, we extensively and analytically studied the system and using the point-splitting regularization skim computed the induced currents (Eq.s (6.20), (6.24), and (6.28)). The net fermion number vanishes due to the C-symmetry. However, the isospin current, J µa T a , and the chiral current, J µ 5 , are non-zero. We found that the dark fermion's backreaction, J = 1 3a δ i a J a i , and its chiral charge, Q 5 = J t 5 , are almost constant during inflation. They are both proportional to the scale of inflation as H 3 , the two sources of CP -violation ξ A and ξ ϕ (corresponding to the VEV of the gauge field and axion), and both increases with the mass of the fermion which explicitly breaks the chiral symmetry in the fermionic sector (Sec. 7.1). In fact, the massless Dirac fermions enjoy conformal and chiral symmetries and have zero (treelevel) induced currents. The massless fermions only have a non-zero chiral charge and therefore backreacts to axion field according to the celebrated Adler-Bell-Jackiw anomaly. Besides, the currents vanish in the absence of the CP -breaking the vacuum. We compared the production of charged scalars with fermions of this setup in Sec. 7.1.1. Interestingly, the vacuum produces more fermions than scalars in two limits; i) when the mass of the corresponding field is higher than Hubble during inflation, i.e., m > H, and ii) in large ξ A limit. Apart from explicit computation for the SU (2)-axion vacuum in Sec. 6, based on symmetries in Sec.s 3 and 7, we explained why the fermion currents should increase with the mass regardless of the source, whether it is the VEV of a SU (2), U (1), or an axion. At first glance, this behavior seems counter-intuitive. However, the present and similar studies which were carried out in de Sitter, are perturbative computations based on assuming a classical background with infinite energy that can always source and support the particle production. However, once the size of backreaction becomes sizable, that approximation is not valid anymore.
We realized that the spontaneous CP -violation during inflation naturally leads to a new and efficient mechanism to generate non-thermal massive dark fermions during inflation. That is regardless of the details of the specific SU (2)-axion model. The size of the fermionic backreaction as well as the present-day density fraction of dark matter put upper bounds on the fermion's mass respectively as 10 2 H (Eq. (7.8)) and (10 10 M Pl H ) 1 4 GeV (Eq. (8.8)). In deriving the latter, we assumed that the heavy fermions are only gravitationally coupled to the rest of the visible and dark sectors. For a GUT scale inflation, the dark fermions can be as heavy as 10 T eV .
The generated dark fermions can be very massive, which depending on the efficiency of their possible feebly interactions with the visible sector, may decay to the visible and lighter dark sectors once the temperature drops below their mass. In other words, these Dirac fermions may replace the right-handed neutrinos in the standard thermal leptogenesis scenario and hence provide a novel natural scenario to explain the observed matter asymmetry in the Universe. However, this possibility, as well as the late time macroscopic properties of the produced (self-interacting) dark fermions as a dark matter sector, are all highly model-dependent. These promising possibilities depend on the value of the gauge coupling, g A , the possible mass of the gauge field after inflation, and the details of the possible feebly interactions of the dark fermions with themselves and the visible sector which is beyond the scope of the current work and we leave for future exhaustive study.

Acknowledgments
Special thanks go to Eiichiro Komatsu for stimulating discussions over the years and bringing the importance of Schwinger effect in the early universe into my attention which led to a series of papers including the current work. I appreciate Kaloian Lozanov for his collaboration in the early stage of the project. I am grateful to David Curtin, Gia Dvali, Elisa Ferreira, Mohammad Khorrami, Yuki Watanabe, and Fabian Schmidt for valuable discussions.

A Mathematical tools
For the sake of self-sufficiency, this appendix reviews mathematical tools that we use throughout this work. This overview includes the definition of the direct sum and product, the spinor covariant derivative, and Whittaker functions, respectively.

A.1 Direct sum and product:
The vector space V, is the direct sum of two subspaces, U 1 and U 2 , as if and only if V = U 1 + U 2 , and U 1 and U 2 are independent. The Kronecker product of two matrices, A m×n and B q×p , is defined as a mp × nq block matrix given by A.2 Spinor covariant derivative: where the spin covariant derivative is with ω µ being the spin-connections and σ αβ = i 4 [γ α , γ β ] being the spinor generators of the Lorentz algebra. The elements of the spin-connection ω αβ µ are given by In FLRW spacetime using the conformal time, the verbeins are and the only non-zero components of the spin connection coefficients are

A.3 Whittaker functions:
Whittaker functions W κ,µ (z) and M κ,µ (z) take the following aymptotic forms in the limit where | arg z |< 3 2 π. Thus, for a complex κ, we have where κ R and κ I are the real and imaginary parts of κ. The Mellin-Barnes integral representation of the Whittaker functions is [63] W κ,µ (z) = e − z 2 2iπ Cα from each other as we show in Fig. 21. In our setup, µ p = i|µ p | is pure imaginary and κ p s is a complex number with a real part equal to |Reκ| = 1 2 .

B Color-spin helicity
This appendix gives a self-contained derivation of the split of our setup as Ψ + ⊕ Ψ − . It is a brief review of appendix B of [50] by the author, which generalized the 2d helicity representation to the 4d color-spin helicity frame (c-helicity).
At this point, it is more insightful to go to the Fourier space and write the system in the Weyl frame. The 8-spinor in the isospin frame can be written in the Weyl frame as whereT 1 is the following 8 × 8 unitary matrix Moreover, the I 2 ⊗ γ α operators transform as where σ α andσ α are and their indices are lowered with the Minkowski metric, i.e. σ α = η αβ σ β . In the Weyl frame, the action (2.13) takes the following form whereL k (τ ) is a 8 × 8 operator given as and Σ 4 is the following 4 × 4 operator As we see in Eq. (B.5), in the absence of µ m , massless case, the system in Eq. (B.4) is decomposed into two independent sub-sectors in terms of the L-and R-handed fields. However, in the massive case, we need to take one step further. The ideal frame would be the one in which Σ 4 is diagonalized. Such frame, if exists, must be made of the common eigenstates of the helicity operator, I 2 ⊗ k i .σ i , and τ i ⊗ σ i . However, these two operators have only two common eigenstates. Thus, Σ 4 is only blockdiagonalizable and reducible up to two subspaces.
As a mathematical tool, we define the color-spin helicity operator for any given 4momentum k α , as in which the first Pauli matrix is a su(2) generator while the second one is the spin operator. The orthonormal eigenstates of this operator are whereǩ α is a four vector given aš Thus,ǩ α is the four momentum of the massless field. The c-helicity h(k α ), and helicity I 2 ⊗k i σ i , have 4 eigenstates in common. However, h(k α ) and τ i ⊗σ i have only two common eigenstates. More precisely, E p s (k α ) with p = ± and s = ± satisfies the eigenstate equations 10) and the orthonormality condition However, only p = + elements of E p s (k α ) are eigenstates of τ i ⊗ σ i as well The E p s (k α ) make an orthonormal basis and hence the following unitary matrix where e p si is the ith element of the E p s (k α ). For each given momentum k α , R k takes the Weyl spinors to their c-helicity frame. Notice that the definition of R k is a bit different than the corresponding matrix in [50]. More precisely, comparing with [50], here we switched the orders of 3rd and 4th columns.
The 8-spinor in the Weyl representation can transformed to the c-helicity state as whereR k is the following 8 × 8 unitary operator In fact,SR k takes the 8-spinor transforms in Weyl frame to the c-helicity framẽ As a result, the Lagrangian operator becomes the following block diagonal 8 × 8 matrix where L ± k, W (τ ) are the following 4 × 4 operations Here,Σ ± (τ, k) are given aš Therefore, the theory in Eq. (B.4) splits into two subsectors as TheL ± k, W (τ ) operators are given as where γ α W s are the gamma matrices in the Weyl representation and notice that Up to now, we split the 8-spinor system into two 4-spinor systems each in a Weyl representation (See Eq. (B.17)). For the quantization purposes, it is more convenient to write the 4-d fields in the Dirac frame. The Weyl 4-spinors can be transformed to their Dirac representation as The gamma matrices in the Weyl representation, γ α W , and Dirac representation, γ α , are related as γ α = Dγ α W D −1 . Therefore, the operatorT 2,k defined as T 2,k ≡DSR k , (B.28) transforms the 8-spinor in the Weyl frame to the final Dirac frame as where now Ψ ± k are each in a Dirac frame. The matrix operators from the isospin to this final extended helicity representation are transformed as Finally, the theory in Eq. (B.4) splits into two subsectors in terms of the plus and minus spinors as This completed the proof that our theory splits into two irreducible representations as (B.37)

C Charge conjugation and Parity
In this appendix, we focus on the discrete symmetries, charge conjugation, and parity. Sec.s C.1 and C.2 investigate the action of the charge conjugation and parity on the theory, respectively.

C.1 Charge conjugation
At this point, we work out the charge conjugation operator and will show it is a symmetry of the theory in Eq. (2.13), i.e.
Notice that the action above is presented in the real space and isospin frame. We define the charge conjugated field asΨ whereC is the charge conjugation operator, and T superscript stands for transpose with respect to Lorentz indices. Under the action of the charge conjugation the fermion charge g A will go to g A,C = −g A which implies Moreover,C is a 8 × 8 matrix which satisfies TheC matrix can be written as Then it is straightforward to check the following equalities Now letsC acts on action in Eq. (C.1). Here for theΨ C field with charge −g A , the action is which can be written in terms ofΨ as where the second to third line has passed by transposition and using the fact that Ψ andΨ anti-commute and in the fourth line we dropped the total derivative terms. By an abuse of notation, it implies That completed our proof that our action is invariant under the charge conjugation and therefore C-symmetric. Equation (C.2) implies that the Fourier modes of the spinor and its charge conjugated field are related asΨ Going to the c-helicity frame, the above reduces to the following form Ψ ± Ck = iγ 2 Ψ ± * −k . (C.11)

(D.2)
Here we recall that µ + = i|µ + |, κ + s = 1 2 + isκ + I andκ + s = − 1 2 + isκ + I . For simplicity, until otherwise stated, we present both κ s andκ s with κ and µ + by µ. Moreover, we drop the symm lim ε→0 while we will keep in mind that all the terms should be symmetrized with respect to the sign of ε, and eventually we take the ε goes to zero limit.
Doing the momentum integral, we arrive at and G(α) is defined as The contours C α and C β can be closed in either the left or right half-plane as far as they separate the poles of Γ( 1 2 + β ± µ) and Γ(−κ − β). (See Fig. 21) In doing the complex integral in G(α) and closing the integral on the right-half plane, we have two types of poles β 1n ≡ −κ * + n, (D.5) where n ∈ N. Therefore, we can write G κ,µ as and G 1 is 16) Doing the integral in Eq. (D.15) and closing the contour on the left-half plane, we have only three poles, α ± = − 1 2 ± µ and α 0 = −1 − κ, which lead to where n ∈ N. The contribution of α 0 n -poles in G 1 is Besides, the contribution of α ± n -poles gives e iκπ /(2π) 2 G 1 | α ± Γ( 1 2 + µ + κ)Γ( (D.20) We computed G 1 and now we turn to compute G 2 .
•• G 2 (τ ; ε) : Now we turn to compute the contribution of the second term inside the square-brackets of G κ,µ in Eq. (D.7), i.e. G 2 as Closing the contour in the right-half plane, the above complex integral, lets call G 2 , can be decomposed as