Axion dark matter in a 3− 3− 1 model J. C. Montero* Universidade Estadual Paulista (UNESP), Instituto de Física Teórica (IFT), São Paulo. R. Dr. Bento Teobaldo Ferraz 271, Barra Funda, São Paulo 01140-070, Brazil Ana R. Romero Castellanos† Instituto de Física Gleb Wataghin—UNICAMP, 13083-859 Campinas, São Paulo, Brazil B. L. Sánchez-Vega‡ Universidade Federal do ABC (UFABC), Centro de Ciências Naturais e Humanas, Av. dos Estados, 5001, 09210-580 Santo André, São Paulo, Brazil and Instituto de Física Gleb Wataghin—UNICAMP, 13083-859 Campinas, São Paulo, Brazil (Received 6 November 2017; published 30 March 2018) Slightly extending a right-handed neutrino version of the 3 − 3 − 1 model, we show that it is not only possible to solve the strong CP problem but also to give the total dark matter abundance reported by the Planck collaboration. Specifically, we consider the possibility of introducing a 3 − 3 − 1 scalar singlet to implement a gravity stable Peccei-Quinn mechanism in this model. Remarkably, for allowed regions of the parameter space, the arising axions with masses ma ≈meV can both make up the total dark matter relic density through nonthermal production mechanisms and be very close to the region to be explored by the IAXO helioscope. DOI: 10.1103/PhysRevD.97.063015 I. INTRODUCTION The impressive observation that almost thirty percent of the energy content of the Universe is due to dark matter (DM) is challenging our understanding of particle physics and cosmology. For a historical review see Ref. [1]. Much effort have been done in order to unravel the nature of DM. Experiments designed to detect weakly interacting massive particles (WIMPs), the, so far, DM candidate paradigm, have failed in providing positive results [2,3]. At the same time, the Large Hadron Collider (LHC) has not been able to produce any signal of a DM candidate, as is the case of the lightest supersymmetric partners of the standard model (SM) neutral particles (gauge or scalar), called neutralinos, or gravitinos (partners of the graviton) [4]. As a consequence of these negative results, it is notice- able the growing interest in studying axions and axionlike particles (ALPs) because they are well motivated alter- natives to WIMPs. Moreover, they can be linked to solutions of still intriguing astrophysical phenomena [5]: (i) ALPs may be the explanation to the TeV photon cosmic transparency if there are gamma ray ↔ ALP oscillations. If so, gamma rays could be converted to ALPs due to the magnetic fields near active galactic nuclei, for instance, traveling “freely” for a long distance to our galaxy and then reconverted into gamma rays in the galactic magnetic fields; (ii) also, ALPs may explain the anomalous energy loss of white dwarfs because from the luminosity of this kind of stars it is inferred that a new energy loss mechanism is needed. In the present scenario, this mechanism could be related to axions or ALPs bremsstrahlung if they directly couple to electrons. All these astrophysical processes constrain the relevant parameters describing axions and ALPs physics. In fact, besides these theoretical arguments for considering axions and/or ALPs, there is also much experimental effort searching for this kind of particles [6]. A variety of experiments have been designed and, in general, they are classified as haloscopes, helioscopes and light-shining through a wall, and most of them are based on the conversion of axions or ALPs into gamma rays in the presence of strong magnetic fields [7]. The axion field was initially introduced as a dynamical solution for the so-called strong CP problem. This problem comes from the extra term which has to be added to the QCD Lagrangian due to the nontrivial structure of the QCD vacuum: *montero@ift.unesp.br †arromero@ifi.unicamp.br ‡bruce.sanchez@ufabc.edu.br Published by the American Physical Society under the terms of the Creative Commons Attribution 4.0 International license. Further distribution of this work must maintain attribution to the author(s) and the published article’s title, journal citation, and DOI. Funded by SCOAP3. PHYSICAL REVIEW D 97, 063015 (2018) 2470-0010=2018=97(6)=063015(12) 063015-1 Published by the American Physical Society https://crossmark.crossref.org/dialog/?doi=10.1103/PhysRevD.97.063015&domain=pdf&date_stamp=2018-03-30 https://doi.org/10.1103/PhysRevD.97.063015 https://doi.org/10.1103/PhysRevD.97.063015 https://doi.org/10.1103/PhysRevD.97.063015 https://doi.org/10.1103/PhysRevD.97.063015 https://creativecommons.org/licenses/by/4.0/ https://creativecommons.org/licenses/by/4.0/ Lθ ¼ θ g2s 32π2 GaμνG̃a μν; whereGaμν is the gluon field strength and G̃a μν its dual. This θ–term violates P, T and CP symmetries and, hence, it induces a neutron electric dipole moment (NEDM). In order to be in agreement with experimental NEDM data the value of the θ parameter must be θ ≲ 0.7 × 10−11 [8]. The strong CP problem is, then, to explain why this parameter is so small. After including weak interactions, the coef- ficient of the GG̃ term changes to θ̄ ¼ θ − arg detMq, where Mq is the quark mass matrix. The Peccei–Quinn (PQ) solution to this problem is implemented by introduc- ing a global U(1) symmetry that must be spontaneously broken and afflicted by a color anomaly. The axion is then the Nambu–Goldstone boson associated to the breaking of that U(1) symmetry, which is now known as the Uð1ÞPQ symmetry. After including the axion field, aðxÞ, the total Lagrangian has a term proportional to the color anomaly NC: LTotal ¼ LSM þ θ̄ g2s 32π2 GaμνG̃a μν þ aðxÞ f̃a=NC g2s 32π2 GaμνG̃a μν þ kineticþ interactions; where f̃a=NC ≡ fa is the axion-decay constant and it is related to the magnitude of the vacuum expectation value (VEV) that breaks the Uð1ÞPQ symmetry. We also have that the divergence of the PQ current, ∂μJ μ PQ, is NC g2s 32π2 GaμνG̃a μν ≠ 0. Hence, the CP violating term GG̃ is now proportional to (θ̄ þ NCaðxÞ=f̃aÞ and it is shown that haðxÞi ¼ −f̃aθ̄=NC minimizes the axion effective potential so that, when the axion field is redefined, aðxÞ → aðxÞ − haðxÞi, the CP violating term GG̃ is no longer present in the Lagrangian, solving in this way the strong CP problem. Although the axion is massless at tree level, it is, in fact, a pseudo-Nambu-Goldstone boson since it gains a mass due to nonperturbative QCD effects related to the Uð1ÞPQ color anomaly. The axion mass and all its couplings are governed by the value of fa. The original conception of the axion was ruled out long ago because fa was thought to be near the electroweak scale, implying in a “visible” axion, in contradiction with laboratory and astrophysical constraints. Few years after the PQ proposal it was realized that for large enough values of fa the axion could be a cold dark matter candidate [9–11]. In fact, for high symmetry breaking scales, the axion is a nonbaryonic extremely weakly-interacting massive particle, stable on cosmological time scales, which makes it a candidate to dark matter. Later in the text we discuss the constraints on fa coming from NEDM, “invisibility” of the axion, and astrophysical data. In order to consider the axion a viable DM candidate we must deal with its relic abundance which strongly depends on the history of the Universe. In particular, the cosmo- logical scenario for the axion production changes signifi- cantly if the PQ symmetry is broken before or after the inflationary expansion of the Universe. The main issue related to the order of these events concerns the axion- production mechanisms. There are production mechanisms due to topological defects, like axionic strings and domain walls, that are comparable to the vacuum misalignment one. Hence, on one hand, if the PQ-symmetry breaking occurs before inflation, inflation will erase these topologi- cal defects. On the other hand, if the PQ-symmetry breaking happens after inflation, it is expected an additional number of axions to be produced due to the decay of the topological defects, affecting directly the relic abundance estimative. In this work we consider axions as DM candidates in the so-called post-inflationary scenario, when the reheating temperature, TR, is high enough to restore the PQ symmetry, TR > TC ∼ fa, which will be broken at a later time, when the temperature of the Universe falls below the critical temperature TC. As we can see, axions present some features with relevant implications not only in particle physics but also in cosmology and it is also a strong indication that physics beyond the SM is in order. In this vein a large variety of models, extensions of the SM, has been proposed. Most of them claim for very appealing achievements relating the DM solution to another yet unsolved issue in particle physics [12–14], as it is the case of the lightness of the active neutrino masses, the smallness of the strong CP violation, or the hierarchy problem, for instance. Among others, a way of introducing new physics is to consider a model with a larger symmetry group. In particular, there is a class of models based on the SUð3ÞC ⊗ SUð3ÞL ⊗ Uð1ÞX gauge group (the so called 3 − 3 − 1 models, for shortness), which are interesting extensions of the SM. In general, these 3 − 3 − 1 models bring welcome features which we review very shortly here. We can take advantage of the larger group representation to choose the matter content in order to introduce new degrees of freedom which are appropriate to implement, for instance, a mecha- nism to generate tiny active neutrino masses, in the lepton sector. The quark sector will also have new degrees of freedom and, depending on the particular representation, the model can have quarks with exotic electric charges or not. The issue of the chiral anomaly cancellation is solved provided we have the same number of triplets and anti- triplets, including color counting. Then, considering that we have the same number of lepton and quark families, say nf, we find that nf must be three or a multiple of three. However, from the QCD asymptotic freedom we find that the number of families must be just three in order to get the correct, negative, sign of the renormalization group β function. Note that, contrarily to the SM, the total number J. C. MONTERO et al. PHYS. REV. D 97, 063015 (2018) 063015-2 of families must be considered altogether in order to get the model anomaly free. Hence, the number of families and the number of colors are related to each other by the anomaly cancellation condition. This fact is a direct consequence of the 3 − 3 − 1 gauge invariance and it can be seen as a hint to the solution to the family replication issue. We can still mention other interesting features: (i) the electric charge quantization does not depend if neutrinos are Majorana or Dirac fermions [15]; (ii) the model described in Refs. [16–18] presents the relation t2 ¼ ðg0=gÞ2 ¼ sin2 θW=ð1 − 4 sin2 θWÞ, which relates the Uð1ÞX and the SUð3ÞL coupling constants, g0 and g, respectively, to the electroweak θW angle. This relations shows a Landau- like pole at some OðTeVÞ, energy scale, μ, for which sin2 θWðμÞ ¼ 1=4 [19], and it would be an explanation to the observed value sin2 θWðMZÞ < 1=4. (iii) The Peccei- Quinn symmetry, usually introduced to solve the strong CP problem, can be introduced in a natural way [20]. In this work we consider a version of a 3 − 3 − 1 model where a gravity stable PQ mechanism can be implemented. We analyze the conditions under which the axion, resulting from the spontaneous breaking of the PQ symmetry in this model, can be considered a dark matter candidate. This work is organized as follows. In Sec. II, we present the general features of the 3 − 3 − 1 model, including its matter content, Yukawa interactions and scalar potential. In Sec. III we show the main steps to make the axion invisible and the PQ mechanism stable against gravitational effects. We also show the axion effective potential from which its mass is derived. In Sec. IV we consider the axion production mechanisms in order to compute its abundance in the Universe. Results for the vacuum misalignment and decay of the string and string-wall system mechanisms are given. In Sec. V we confront the predictions from the previous section with the observational constraints, coming mainly from the Planck-collaboration results for the DM abundance, the NEDM data and direct axion searches, in order to constrain the parameter space of the model. Section VI is devoted to our final discussions and conclusions. II. BRIEFLY REVIEWING THE MODEL We consider the 3 − 3 − 1 model with right-handed neutrinos, Na, in the same multiplet as the SM leptons, νa and ea. In other words, in this model all of the left- handed leptons, FaL ¼ ðνa; ea; Nc aÞTL with a ¼ 1, 2, 3, belong to the same ð1; 3;−1=3Þ representation, where the numbers inside the parenthesis denote the quantum numbers of SUð3ÞC, SUð3ÞL and Uð1ÞX gauge groups, respectively. This model was proposed in Refs. [21,22] and it has been subsequently considered in Refs. [20,23–31]. It shares appealing features with other versions of 3 − 3 − 1 models [16–18,32–36]. Furthermore, the existence of right- handed neutrinos allows mass terms at tree level, but it is necessary to go to the one-loop level to obtain neutrino masses in agreement with experiments [26]. The remaining left-handed fermionic fields of the model belong to the following representations Quarks∶ QL ¼ ðu1; d1; u4ÞTL ∼ ð3; 3; 1=3Þ; ð1Þ QbL ¼ ðdb; ub; dbþ2ÞTL ∼ ð3; 3̄; 0Þ; ð2Þ where b ¼ 2, 3; and “∼” means the transformation proper- ties under the local symmetry group. Additionally, in the right-handed field sector we have Leptons∶ eaR ∼ ð1; 1;−1Þ; ð3Þ Quarks∶ usR ∼ ð3; 1; 2=3Þ; dtR ∼ ð3; 1;−1=3Þ; ð4Þ where a ¼ 1, 2, 3; s ¼ 1;…; 4 and t ¼ 1;…; 5. In order to generate the fermion and boson masses, the SUð3ÞC ⊗ SUð3ÞL ⊗ Uð1ÞX symmetry must be spontane- ously broken to the electromagnetic group, i.e., to the Uð1ÞQ symmetry, whereQ is the electric charge. To do this, it is necessary to introduce, at least, three SUð3ÞL triplets, η, ρ, χ, as shown in Ref. [30], which are given by η ¼ ðη01; η−2 ; η03ÞT ∼ ð1; 3;−1=3Þ; ρ ¼ ðρþ1 ; ρ02; ρþ3 ÞT ∼ ð1; 3; 2=3Þ; ð5Þ χ ¼ ðχ01; χ−2 ; χ03ÞT ∼ ð1; 3;−1=3Þ: ð6Þ Once these fermionic and bosonic fields are introduced in the model, we can write the most general Yukawa Lagrangian, invariant under the local gauge group, as follows LYuk ¼ Lρ Yuk þ Lη Yuk þ Lχ Yuk; ð7Þ with Lρ Yuk ¼ αtQ̄LdtRρþ αbsQ̄bLusRρ� þ Yaa0ϵijkðF̄aLÞiðFa0LÞcjðρ�Þk þ Y0 aa0F̄aLea0Rρ þ H:c:; ð8Þ Lη Yuk ¼ βsQ̄LusRηþ βbtQ̄bLdtRη� þ H:c:; ð9Þ Lχ Yuk ¼ γsQ̄LusRχ þ γbtQ̄bLdtRχ� þ H:c:; ð10Þ where ϵijk is the Levi-Civita symbol and a0, i; j; k ¼ 1, 2, 3 and a, b, s, t are in the same range as in Eq. (3). It is also straightforward to write down the most general scalar potential consistent with gauge invariance and renormalizability as AXION DARK MATTER IN A 3 − 3 − 1 MODEL PHYS. REV. D 97, 063015 (2018) 063015-3 Vðη; ρ; χÞ ¼ VZ2 ðη; ρ; χÞ þ VZ2 ðη; ρ; χÞ; ð11Þ with VZ2 ðη; ρ; χÞ ¼ −μ21η†η − μ22ρ †ρ − μ23χ †χ þ λ1ðη†ηÞ2 þ λ2ðρ†ρÞ2 þ λ3ðχ†χÞ2 þ λ4ðχ†χÞðη†ηÞ þ λ5ðχ†χÞðρ†ρÞ þ λ6ðη†ηÞðρ†ρÞ þ λ7ðχ†ηÞðη†χÞ þ λ8ðχ†ρÞðρ†χÞ þ λ9ðη†ρÞðρ†ηÞ þ ½λ10ðχ†ηÞ2 þ H:c:�; ð12Þ VZ2 ðη; ρ; χÞ ¼ −μ24χ†ηþ λ11ðχ†ηÞðη†ηÞ þ λ12ðχ†ηÞðχ†χÞ þ λ13ðχ†ηÞðρ†ρÞ þ λ14ðχ†ρÞðρ†ηÞ þ λ15ϵijkηiρjχk þ H:c: ð13Þ We have divided the total scalar potential Vðη; ρ; χÞ in two pieces, VZ2 ðη; ρ; χÞ, invariant under the Z2 discrete sym- metry (χ → −χ, u4R → −u4R, dð4;5ÞR → −dð4;5ÞR, and all the other fields even by the symmetry), and VZ2 ðη; ρ; χÞ, which breaks Z2. This discrete symmetry is motivated by the implementation of the PQ mechanism as shown below. It is well known that the minimal vacuum structure needed to give masses to all the particles in the model is hρi ¼ 1ffiffiffi 2 p ð0; vρ0 2 ; 0ÞT; hηi ¼ 1ffiffiffi 2 p ðvη0 1 ; 0; 0ÞT; hχi ¼ 1ffiffiffi 2 p ð0; 0; vχ0 3 ÞT; ð14Þ which correctly reduces the SUð3ÞC ⊗ SUð3ÞL ⊗ Uð1ÞX symmetry to the Uð1ÞQ one. In principle, the remaining neutral scalars, η03 and χ 0 1, can also gain VEVs. However, in this case, dangerous Nambu-Goldstone bosons can arise in the physical spectrum, as shown in Ref. [37]. In this paper, we are going to consider only the minimal vacuum structure given in Eq. (14). III. IMPLEMENTING A GRAVITY STABLE PQ MECHANISM The key ingredient to implement the PQ mechanism is the invariance of the entire Lagrangian under a global Uð1Þ symmetry, called Uð1ÞPQ, which must be both afflicted by a color anomaly and spontaneously broken [38–41]. In general, the implementation of the PQ mecha- nism in the 3 − 3 − 1 models is relatively straightforward [20,28]. In particular, in Ref. [20] a gravitationally stable PQ mechanism for the model considered here is success- fully implemented. We are going to review its main results for completeness. First of all, we search for all Uð1Þ symmetries of the Lagrangian given in Eqs. (7) and (11). Doing so, we find only two symmetries, Uð1ÞX and Uð1ÞB, which clearly do not satisfy the two minimal conditions required for the Uð1ÞPQ symmetry. See Table I for the quantum number assignments of the fields for these symmetries. In other words, the Uð1ÞPQ is not naturally allowed by the gauge symmetry. However, if the Lagrangian is slightly modified by imposing a Z2 discrete symmetry such that χ → −χ, u4R → −u4R, dð4;5ÞR → −dð4;5ÞR, all terms in VZ2 ðη; ρ; χÞ are forbidden. In addition, the Yukawa Lagrangian inter- actions given in Eqs. (8)–(10) are slightly modified to Lρ Yuk ¼ αaQ̄LdaRρþ αbaQ̄bLuaRρ� þ Yaa0εijkðF̄aLÞiðFbLÞcjðρ�Þk þ Y0 aa0F̄aLea0Rρ þ H:c:; ð15Þ Lη Yuk ¼ βaQ̄LuaRηþ βbaQ̄bLdaRη� þ H:c:; ð16Þ Lχ Yuk ¼ γ4Q̄Lu4Rχ þ γbðbþ2ÞQ̄bLdðbþ2ÞRχ� þ H:c: ð17Þ Consequently, with the imposition of this Z2 symmetry a Uð1ÞPQ symmetry is automatically introduced with the charges given in Table II. As η, ρ, χ get VEVs, an axion appears in the physical spectrum. However, it is a visible axion because the Uð1ÞPQ symmetry is actually broken by vρ0 2 , which is upper bounded by the value of vSM ≃ 246 GeV, as shown in Refs. [20,37]. Hence, this scenario is ruled out [42]. Nevertheless, a singlet scalar, ϕ ∼ ð1; 1; 0Þ, can be intro- duced in order to make the axion invisible. Its role is to break the PQ symmetry at an energy scale much larger than the electroweak one. This field does not couple directly to quarks and leptons, however it couples to the scalar triplets, η, ρ and χ, through Hermitian terms and the non-Hermitian term λPQϵ ijkηiρjχkϕ, from which it gets a PQ charge equal to 6, cf. Table II. Notice that this term is allowed as long as the ϕ field is odd under theZ2 symmetry, i.e.,Z2ðϕÞ ¼ −ϕ. TABLE I. The U(1) symmetries of the Lagrangian given by Eqs. (7) and (11). QL QiL (uaR, u4R) (daR, dð4;5ÞR) FaL eaR ρ (χ, η) Uð1ÞX 1=3 0 2=3 −1=3 −1=3 −1 2=3 −1=3 Uð1ÞB 1=3 1=3 1=3 1=3 0 0 0 0 TABLE II. The Uð1ÞPQ charges in the model with a Z2 discrete symmetry such that χ→−χ, u4R → −u4R, and dð4;5ÞR → −dð4;5ÞR. QL QiL (uaR, u4R) (daR, dð4;5ÞR) FaL eaR ρ (χ, η) Uð1ÞPQ −2 2 0 0 1 3 −2 −2 J. C. MONTERO et al. PHYS. REV. D 97, 063015 (2018) 063015-4 Although the Z2 discrete symmetry apparently introdu- ces the PQ mechanism in the model, there are two issues with it. First, the Z2 and gauge symmetries allow some renormalizable terms in the scalar potential, such as ϕ2, ϕ4, ρ†ρϕ2, η†ηϕ2, χ†χϕ2, that explicitly violate the PQ sym- metry in an order low enough to make the PQ mechanism ineffective. Second, since the PQ symmetry is global, it is expected to be broken by gravitational effects [43,44]. Thus, a mechanism to stabilize the axion solution has to be introduced. As usual, the entire Lagrangian is considered to be invariant under a ZD discrete gauge symmetry (anomaly free) [20,28,45–48] and, in addition, this symmetry is supposed to induce the Uð1ÞPQ symmetry. For ZD≥10 it is found that all effective operators of the form ϕN=MN−4 Pl (where N ≥ D is a positive integer and MPl is the reduced Planck mass) that can jeopardize the PQ mechanism are suppressed. In particular, in Ref. [20] two different sym- metries, Z10 and Z11, were found to stabilize the PQ mechanism for the Lagrangian given by Eqs. (11), (15)– (17). The specific charge assignments for these symmetries are shown in Table III. Note that the term λ15ϵijkηiρjχk in the scalar potential is prohibited by both of these discrete symmetries and it must be removed from the entire Lagrangian. We remark that both theZ10 andZ11 discrete symmetries in Table III are anomaly free. This type of discrete symmetry is known as gauge discrete ZN symmetry and it is assumed to be a remnant of a gauge (local) symmetry valid at very high energies, [45]. The anomaly-free con- ditions are necessary in order to truly protect the PQ mechanism against gravity effects [46,49–51], Specifically, these discrete symmetries satisfy A3CðZNÞ ¼ A3LðZNÞ ¼ 0 Mod N=2, where A3C and A3L are the ½SUð3ÞC�2 × ZN , ½SUð3ÞL�2 × ZN anomalies, respectively. Other anomalies, such as Z3 N , do not give useful low energy constraints because these depend on some arbitrary choices concerning to the full theory. In particular, the Z3 N anomaly depends on the fermions which get masses at very high energy and are integrated out in the low-energy Lagrangian. All the details of these anomaly conditions applied to the 3 − 3 − 1model can be found in Ref. [20]. In both cases, the axion, aðxÞ, is the phase of the ϕ field, i.e., ϕðxÞ ∝ exp ðiaðxÞ=f̃aÞ, which implies f̃a ≈ vϕ. As it is well known, to make the axion compatible with astro- physical and cosmological considerations, the axion-decay constant fa (related to f̃a by fa ¼ f̃a=NC ¼ f̃a=NDW, with NDW being the number of domain walls in the theory. In this model we have NC ¼ NDW ¼ 3), must be in the range 109 GeV≲ fa ≲ 1012 GeV (we are assuming a post- inflationary PQ symmetry breaking scenario). Note that this high value of fa ¼ f̃a=NC ≈ vϕ=NC ≫ vρ0 2 ; vη0 1 ; vχ0 3 , justifies the approximation in the form of axion eigenstate. It is also important to remember that in this model v2 ρ0 2 þ v2 η0 1 ¼ v2SM and vχ0 3 is expected to be at the TeVenergy scale. Now, we can go further calculating the axion mass, ma. In this model, the axion gains mass because the Uð1ÞPQ symmetry is both anomalous under the SUð3ÞC group and explicitly broken by gravity-induced operators, gϕN=MN−4 Pl (with g ¼ jgj exp iδ). These operators have a high dimen- sion (N ≥ 10) because of the protecting Z10 or Z11 discrete symmetries, as shown in Table III. These two effects induce an effective potential for the axion, Veff , from which it is possible to determine the axion mass. In more detail, as the Uð1ÞPQ symmetry is anomalous, we will have a VPQ term in the effective potential, which can be written as VPQ ¼ −m2 πf2π � 1 − 4mumd ðmu þmdÞ2 sin2 � aðxÞ 2fa �� 1=2 ; ð18Þ where mπ ≃ 135 MeV and fπ ≃ 92 MeV are the mass and decay constant of the neutral pion, respectively; mu and md are the masses of the up and down quarks. Note that VPQ has a minimum when haðxÞi=fa ¼ 0, which solves the strong CP problem in the usual way. However, because of the PQ symmetry is also explicitly broken by gravity effects, the effective potential gets another term, Vgravity, which reads Vgravity ≃ − jgjvNϕ 2N=2−1MN−4 Pl cos � NaðxÞ f̃a þ δD � ; ð19Þ where N ¼ 10, 11 for Z10 and Z11, respectively. The phase δD inside the trigonometric function can be written as δD ¼ δ − Nθ̄; ð20Þ where δ is the phase of the g coupling constant and θ̄ is the parameter which couples to the gluonic field strength and its dual. This extra term in the scalar potential, Eq. (19), has two important consequences. First, it induces a shift in the value of haðxÞifa where Veff has a minimum. Expanding Veff ¼ VPQ þ Vgravity in powers of haðxÞi fa , we find that in the minimum, the axion VEV satisfies jhaðxÞij fa ���� min ≃ NjgjNN−1 DW 2 N 2 −1 ð fa MPl ÞN−2M2 Pl sin δD m2 πf2π f2a mumd ðmuþmdÞ2 þ N2jgjNN−2 DW 2 N 2 −1 ð fa MPl ÞN−2M2 Pl cos δD ; ð21Þ TABLE III. The charge assignment for ZD that stabilizes the PQ mechanism in the considered 3 − 3 − 1 model. QL QiL (uaR, u4R) (daR, dð4;5ÞR) FaL eaR ρ (χ, η) ϕ Z10 þ7 þ5 þ1 þ1 þ7 þ1 þ6 þ6 þ2 Z11 þ7 þ6 þ1 þ1 þ8 þ2 þ6 þ6 þ4 AXION DARK MATTER IN A 3 − 3 − 1 MODEL PHYS. REV. D 97, 063015 (2018) 063015-5 where we have used vϕ ≈ f̃a ¼ NDWfa. Note that for jgj ¼ 0 (or for δD ¼ 0) we have that haðxÞi fa ¼ 0 in the minimum, as it should be to solve the strong CP problem. However, in the general case, the value of haðxÞi fa does not satisfy the NEDM constraint [8], which imposes haðxÞi fa ¼ θ̄ ≲ 0.7 × 10−11: ð22Þ In addition, Vgravity brings a mass contribution for the axion, ma;gravity. From Eq. (19) we obtain m2 a;gravity ¼ N2jgjNN−2 DW 2 N 2 −1 � fa MPl � N−2 M2 Pl cos δD: ð23Þ This contribution can, in general, be much larger than the well-known axion-mass term coming from the QCD non- perturbative terms, Eq. (18), m2 a;QCD ¼ m2 πf2π f2a mumd ðmu þmdÞ2 : ð24Þ Thus, in order to maintain the axion mass stable, we are going to look for values of the parameters jgj, fa and δD for N ¼ 10, 11 that both satisfy the NEDM constraint and leave the axion mass stable ðma;QCD ≳ma;gravityÞ. Before closing this section, it is important to remark that although the 3 − 3 − 1 model considered in this paper has additional contributions to CP-violating processes that in principle can contribute to the NEDM, these do not require tuning the model parameters at the same order of the θ̄ parameter as it was correctly estimated in Ref. [20]. Roughly speaking, the dominant contribution to the up-quark electric dipole moment, deu, coming from the interchange of the χ scalar is of order deujmu≪mu4 ;mχ ≈ ejγ4·γbðbþ2Þj sin α 48π2 mu4 m2 χ KðrÞ, where sin α is the sine of the CP-violating phase, α, and KðrÞ ¼ 1 2r − 1 r2 þ 1 r3 ln ð1þ rÞ, with r ¼ m2 u4 m2 χ − 1, and where mu is the up-quark mass; mu4 andmχ are the exotic quark and scalar masses, respectively. For reasonable Yukawa couplings (γ4, γbðbþ2Þ) and CP-violating phases, and for mu4 and mχ masses of order of TeV, the den ∼ 4 3 ded − 1 3 deu ≈OðdeuÞ is in agreement with experiments without requiring a strong fine-tuning of the parameter of the model [20]. IV. REVIEWING THE NONTHERMAL PRODUCTION OF AXION DARK MATTER For the postinflationary fa values considered here, cold dark matter in the form of axions can be produced by three different processes: the misalignment mechanism [52], where the axion field oscillates about the minimum of its potential, trying to decrease the energy after the breaking of the PQ symmetry; and the decay of one-dimensional (global strings [53]) and two-dimensional (domain walls [54]) topological defects, which appear after breaking this symmetry. Now, we will briefly review the general expres- sions for the axion relic density in these three mechanisms following Ref. [55]. A. Misalignment mechanism The equation of motion for the axion field a in a homogeneous and isotropic Universe, is of the type of a damped harmonic oscillator with a natural frequency equal to the axion mass. In this case, taking into account nonperturbative effects of QCD at finite temperature and considering the interacting instanton liquid model (IILM) [56], the axion mass depends on the temperature as [57] m2 aðTÞ ¼ cT Λ4 QCD f2a � T ΛQCD � −n ; ð25Þ where the values of the parameters are cT ¼ 1.68 × 10−7, n ¼ 6.68 and ΛQCD¼ 400 MeV [57]. This dependence, is valid in the regimewhere the axion mass at temperature T is less than its value at temperature zero, given by mað0Þ2 ¼ c0 Λ4 QCD f2a , where c0 ¼ 1.46 × 10−3, which leads to a mini- mum temperature ∼103 MeV for the validity of the fit. The temperature Tosc at which the axion field begins to oscillate is given by [55] Tosc ¼ 2.29 GeV � g�ðToscÞ 80 � − 1 4þn � fa 1010 GeV � − 2 4þn × � ΛQCD 400 MeV � ; ð26Þ where g�ðToscÞ is the number of relativistic degrees of freedom at temperature Tosc. Eq. (26) is valid for temper- atures greater than 103MeV, where Eq. (25) holds, and it is also assumed a not too strong dependence on the temper- ature of g�, which, for the range 109 GeV < fa < 1012 GeV analyzed in this work, varies between 80 and 85 [58], what would change the abundance of axion dark matter by a factor of ≈1.02. Once the adiabatic condition is satisfied, both the entropy and the number of axions with momentum zero per comoving volume are con- served [9], and it is possible to obtain the dark matter abundance [55] Ωa;mish2 ¼ 4.63 × 10−3 � fa 1010 GeV �6þn 4þn ; ð27Þ J. C. MONTERO et al. PHYS. REV. D 97, 063015 (2018) 063015-6 where g�ðToscÞ ¼ 80 and ΛQCD¼ 400 MeV have been used. B. Decay of global strings Global strings are the first of the topological defects that appear after the breaking of the Uð1ÞPQ symmetry at T ≲ vϕ because the field ϕ (with PQ charge equal to 6 in the 3 − 3 − 1 model considered here) acquires a VEV jhϕij ¼ vϕ [55,59]. Actually, the breaking of the PQ symmetry leads to the formation of a densely knotted network of cosmic axion strings, which oscillate under their own tension, losing their energy by radiating axions [60]. The radiation process lasts from the PQ-symmetry breaking time to the QCD phase transition time. Using results of numerical studies which provide the time dependence of ρstring (energy density of strings) and ρa;string (energy density of axions produced by the string decays), it is possible to obtain the nowadays abundance of radiated axions [61,62], Ωa;stringh2 ¼ αN2 DW × � fa 1010 GeV �6þn 4þn ; ð28Þ with α ¼ ð7.3� 3.9Þ × 10−3, g�ðToscÞ ¼ 80 and ΛQCD¼ 400 MeV. NDW ¼ 3 is the number of domain walls in this model, and n ¼ 6.68 is the same parameter that appears in Eq. (25). C. Decay of string-wall systems In the 3 − 3 − 1 model considered, a Z3 subgroup remains after the breaking of the Uð1ÞPQ symmetry, which makes the vacuum manifold to be made of several disconnected components. When the temperature of the Universe lies between the electroweak and QCD phase transition energy scales, domain walls appear as a conse- quence of breaking this Z3 discrete symmetry. These domain walls are attached by strings and occur at the boundaries between regions of space-time where the value of the field ϕ is different. These inhomogeneities of space-time are in tension with the assumptions of standard cosmology. So, it is necessary that these domain walls decay at a certain time after being formed [63]. Actually, the domain walls bounded by strings begin to oscillate and eventually, when their tensions are greater than the tensions of the strings, their annihilations lead to axion production [64,65]. The energy density of domain walls can overclose the Universe due to its dependence on the inverse of the square of the scale factor, R, which decreases at a slower rate than the corresponding to matter, ρ∼R−3, and radiation, ρ∼R−4. In our case, this problem is solved by the introduction of a Planck-suppressed operator in the effective potential for the axion field a, parametrized as in Eq. (19). The current axion abundance is given by the expression [55,66]: Ωa;wallh2 ¼ 1.23 × 10−6½7.22 × 103� 32pβ � 2p − 1 3 − 2p � × � N4 DW � 1 − cos 2πN NDW �� 1− 3 2p × jgj1− 3 2p � Ξ 10−52 � 1− 3 2p � fa 1010 GeV � 4þ3ð4p−16−3nÞ 2pð4þnÞ ; ð29Þ where Ξ ¼ 1 2 N 2 ð vϕ MPl ÞN−4, and β ¼ 1.65� 0.47 is a parameter obtained from numerical simulations. Finally, we will refer to the case p ¼ 1 as the exact scaling, and p ≠ 1 as the deviation from scaling. From here on, we use p ¼ 0.926 for the deviation from scaling case, since it is the suggested value by numerical simulations [55]. In order to conclude this section, we have seen that axions can be produced by three different non-thermal mechanisms, which leads to the result that the total abun- dance of axions in the Universe can be written as the sum of all these contributions, Eqs. (27), (28) and (29), i.e., Ωah2 ¼ Ωa;mish2 þΩa;stringh2 þ Ωa;wallh2: ð30Þ The total dark matter abundance due to axions is upper bounded by the observational constraint on the current relic density ΩPlanck DM h2 ¼ 0.1197� 0.0066 (at 3σ) as reported by the Planck Collaboration [67]. In the next section, we will analyze the behavior of each contribution to the total abundance, in order to establish a suitable region of parameters for the model analyzed in this work. V. CONSTRAINING THE NONTHERMAL PRODUCTION OF AXION DARK MATTER In general, the total dark matter relic density due to axions in this 3 − 3 − 1 model depends on fa, g; NDW and ZN . The dependence on fa, g, and NDW is direct because Ωa;mis;Ωa;string and Ωa;wall explicitly depend on these parameters. Nevertheless, the dependence on ZN is indirect. Roughly speaking, this discrete symmetry con- strains the order of the dominant gravity-induced operator gϕN=MN−4 Pl . In other words, the discrete symmetry sets the exponent N which directly affects the total dark matter due to axions. Actually, we have two discrete symmetries, Z10 and Z11 (see Table III), that stabilize the PQ mechanism, which implies that there are two cases to be considered, N ¼ 10 and N ¼ 11. On the other hand, the domain wall parameter, NDW, is set to be equal to 3 by the PQ symmetry and the matter content in the model. Thus, we are interested in knowing if the model with Z10 or/and Z11 symmetry provides the total dark matter reported by the Planck collaboration [67] when fa, g, take their allowed values, without conflicting with the constraints on the axion phenomenology. AXION DARK MATTER IN A 3 − 3 − 1 MODEL PHYS. REV. D 97, 063015 (2018) 063015-7 In order to do that, it is convenient, first, to study separately the behavior of the three axion production mechanisms which results are shown in Fig. 1. Specifically, the cyan and black lines show the axion abundances produced by misalignment and global string decay mechanisms, respectively. On the other hand, the blue lines show the abundance of axion dark matter due to the decay of domain wall systems for N ¼ 10 and N ¼ 11, calculated for the coupling constant value jgj ¼ 1. Two shaded regions are also shown: the light red one corre- sponds to the exclusion region coming from the constraint of the over closure of the Universe [67], and the yellow region gives the possible interval for the axion decay constant fa, for which no over abundance of axions from decay of global strings or domain walls is produced. Finally, the dark green line corresponds to the total abundance of axions, Ωah2, as given by Eq. (30), obtained for the caseN ¼ 10 and jgj ¼ 1. The case forN ¼ 11 is not shown because for all the considered values of fa the axion relic density is overabundant. From Fig. 1 some conclusions are straightforward. First, Ωa;mis and Ωa;string grow when fa grows. Thus, in principle, these are dominant for the greater values of fa (5.3 × 109 GeV≲ fa ≲ 1.7 × 1010 GeV). However, the misalignment mechanism is always subdominant because Ωa;string has an extra N2 DW ¼ 9 global factor. Indeed, the misalignment mechanism contributes at most by ≈7% for the total dark matter density. In contrast, Ωa;wall is decreasing with fa and thus it dominates Ωa for the smaller values of fa (3.6 × 109 GeV≲ fa ≲ 5.3 × 109 GeV). That can be understood realizing that the domain-wall time decay is larger for smaller fa values, making the domain wall more stable and, in this way, explaining why this mechanism contributes more for the axion relic density when fa is smaller. The opposite behavior of Ωa;string and Ωa;wall allow to set an upper and lower bound on fa. For jgj ¼ 1, fa is constrained to be 3.6×109 GeV ma;gravity [see the discussion near Eq. (23)]. Becausema;gravity is directly proportional to jgj and fN−2 a , cf. Eq. (23), and m2 a;QCD is inversely proportional to f2a, cf. Eq. (24), the forbidden region, denoted by the light red color, is in the top right part of the fa − jgj plane. In addition, in Fig. 2 are shown three dark red lines which correspond to the NEDM constraint, given by Eq. (22), for different values of δD. It is important to realize that δD values of order one (not shown) do not give allowed regions in the parameter space. It is necessary to allow δD ≲ 10−5 in order to have nonexcluded regions which are below the lines. In particular, we calculate the maximum values of δD that give allowed regions in the parameter space. The corresponding results, in the cases of exact scaling (p ¼ 1) and deviation from scaling (p ¼ 0.926), are δD ¼ � ð0.4–4.1Þ × 10−5 Exact scaling; ð2.9–9.5Þ × 10−6 Deviation from scaling: ð31Þ FIG. 1. Relic density of nonthermal axion dark matter in the 3 − 3 − 1model, assuming exact scaling, p ¼ 1, and jgj ¼ 1. The central values of the parameters in Eqs. (28) and (29) together with NDW ¼ 3 have been used. The vertical dashed lines limit regions with over production of axions by decay of domain walls (left line) and strings (right line), while the horizontal red line is the experimental constraint Ωah2 ¼ ΩPlanck DM h2. J. C. MONTERO et al. PHYS. REV. D 97, 063015 (2018) 063015-8 These values are obtained by taking jgj ¼ ffiffiffiffiffiffi 4π p , and con- sidering the uncertainties in the parameters of the three axion production mechanisms. Lower values of jgj would require higher tuning on the δD parameter, with values of the order 10−8 as shown in Fig. 2. In general, for jgj fixed, the tuning on δD depends on the decay constant fa and the mechanism of axion darkmatter production: if the decay of domain walls was dominant (left side of the curves), the tuning would be less severe than if the production by string decay (right side of the curves) was the dominant one. Also, in Fig. 2 is shown that for a δD small enough in order to satisfy the NEDM condition, and for a given jgj value between 5 × 10−2 and ffiffiffiffiffiffi 4π p , there are two separated regions for fa where axions can make up the total DM relic density. For instance, taking jgj ¼ ffiffiffiffiffiffi 4π p and considering the uncertainties in the parameters, these regions and their corresponding axion masses for the exact scaling case, are fa ≈ � ð2.8–3.5Þ× 109 GeV→ma ≈ ð1.7–2.1Þ× 10−3 eV ð1.1–1.2Þ× 1010 GeV→ma ≈ ð5–5.4Þ× 10−4 eV ð32Þ In the first range for fa the production of dark matter is mainly through the decay of domain walls, while in the second range it is due to the decay of strings. Taking smaller values for jgj, will lead to more stringent intervals for both fa and ma. For the case of deviation from scaling, we find fa ≈ ð3.4–3.6Þ × 109 GeV, corresponding to ma ≈ ð1.7–1.8Þ × 10−3 eV, when the domain walls decay is the leading production mechanism, and fa ≈ ð1.1–1.2Þ× 1010 GeV, leading to ma ≈ ð5–5.4Þ × 10−4 eV, for the string decay as the dominant contribution. Finally, for values of jgj of order one, we can make predictions regarding the observability of axion in current and/or future experiments. Specifically, the axion coupling to two photons, gaγγ , depends on the fa decay constant, the electromagnetic and color anomaly coefficients, E and NC, respectively. It is known that these anomaly coefficients are completely determined by the fermion content and the (a) (b) FIG. 2. Observational constraints on the parameter space fa − jgj in the 3 − 3 − 1 model, assuming exact scaling (a) and deviation from scaling (b). These plots correspond to the Z10 discrete symmetry, and NDW ¼ 3. The shaded regions in light red and light blue correspond to regions of the parameter space where the constraints given by ma;QCD > ma;gravity and Ωah2 ≤ ΩPlanck DM h2 are violated, respectively. Moreover, the regions above the straight red lines correspond to the exclusion regions set by the NEDM condition, as given by Eqs. (21) and (22), for three different choices of the δD parameter. FIG. 3. Projected sensitivities of different experiments in the search for axion dark matter. The green regions show sensitivities of light-shining-through-wall experiments like ALPS-II [70], of the helioscope IAXO [69], of the haloscopes ADMX and ADMX-HF [71,72]. The yellow band corresponds to the generic prediction for axion models in QCD. In addition, the two (one) thick red (blue) lines stand for the predicted mass ranges and coupling to photons in this model, for jgj ¼ 0.1 (jgj ¼ 1), where axions make up the total DM relic density. AXION DARK MATTER IN A 3 − 3 − 1 MODEL PHYS. REV. D 97, 063015 (2018) 063015-9 Uð1ÞPQ charges of the model, cf. Table (II). Standard calculations for anomaly coefficients [41,68] furnish E ¼ −4 and NC ¼ 3. With this information, we can go further plotting, in Fig. 3, gaγγ as a function of ma for the regions where axions make up the total dark matter relic density and for two different values of jgj, specifically jgj ¼ 0.1 and jgj ¼ 1. This figure clearly shows two allowed regions for jgj ¼ 0.1: ma ≈ ð0.4–0.6Þ × 10−3 eV with gaγγ ≈ ð4.5–5.9Þ × 10−13 GeV−1 and ma ≈ ð0.9–1.3Þ × 10−3 eV with gaγγ ≈ ð1.1–1.6Þ × 10−12 GeV−1, and one region for jgj ¼ 1: ma ≈ ð1.4–1.8Þ × 10−3 eV with gaγγ ≈ ð1.8–2.2Þ × 10−12 GeV−1. The reason why there is only one region for larger jgj values is that the gravitational mass grows with jgj and thus, it conflicts with the condition ma;QCD ≫ ma;gravity for lower axion masses. Moreover, it is notable that for the range with larger masses (blue line), the axion parameters of this 3 − 3 − 1 model are very close to the projected region which is going to be explored by the IAXO experiment [66,69]. VI. CONCLUSIONS In this work, we consider a version of an alternative electroweak model based on the SUð3ÞL ⊗ Uð1ÞX gauge symmetry, the so called 3 − 3 − 1 models, when the color gauge group is added. For this version, which includes right-handed neutrinos, it is shown in Ref. [20] that the PQ mechanism for the solution of the strong CP problem can be implemented. In this implementation, the axion, the pseudo Nambu-Goldstone boson that emerges from the PQ-symmetry breaking, is made invisible by the introduc- tion of the scalar singlet ϕ ∼ ð1; 1; 0Þ whose VEV, vϕ ≈ f̃a, is much larger than vSM, and any other VEV in the model. Moreover, the axion is also protected against gravitational effects, that could destabilize its mass, by a discrete ZN symmetry, with N ¼ 10, 11. Once we have set this consistent scenario, we investigate the capabilities of this axion, produced in the framework of this particular 3 − 3 − 1 model, to be a postinflationary cold dark matter candidate. We started focusing in the axion-production mechanisms. As it was explained in the previous section, from Fig. 1 we see that the vacuum misalignment mechanism does not dominate the DM relic abundance, and, if it was the only production mechanism in action, an upper bound for fa could be set by imposing that it should account for all the DM abundance, i.e., Ωa;mish2 ¼ ΩPlanck DM h2, and we would find the corresponding value fa ≈ 1.5 × 1011 GeV, for the parameters determined by the model, in this caseNDW ¼ 3. However, there are two other more efficient mechanisms due to the decay of topological defects: cosmic strings and domain walls. As the curves for Ωa;stringh2 and Ωa;wallh2 grow in opposite directions, relatively to the fa values, we can determine an upper bound and a lower bound for fa by imposing the total Ωah2 matches the observed Planck results. This is the case when we add up all the contributions for N ¼ 10, and we find 3.6 × 109 GeV < fa < 1.7 × 1010 GeV. However, we would like to stress that this is not the case for N ¼ 11. For N ¼ 11 there is no value of fa for which the addition of the partial abundances lies below the observed result. It means that theZ11, which possesses the good quality of stabilizing the axion, is not appropriate for the axion-production issue since it makes the domain wall mechanism too efficient and overpopulates the Universe. As it can be seen from Fig. 1, for any fixed allowed value of jgj, there are two values of fa that are in agreement with the value of ΩPlanck DM h2. In fact they are regions, if we take into account the uncertainties following the discussion in the previous section for Fig. 2. Outside these regions, the axion abundance will be a fraction of ΩPlanck DM h2. See the solid dark green curve in Fig. 1 for jgj ¼ 1. If this happens to be the case, i.e., if these predicted regions are somehow excluded, by future experimental data for the axion mass value, for instance, then, another kind of DM will be needed. We have also found special values for δD, ð0.4–4.1Þ × 10−5, by requiring the minimal compatible intersection region between the curves that obey the NEDM and ΩPlanck DM h2 constraints. This value was obtained considering the maximum value of jgj, i.e., jgj ¼ ffiffiffiffiffiffi 4π p , cf. Fig. 2(a). For lower values of jgj, higher tuning on δD is required. However, it seems unnatural to require severe levels of tuning on δD, since for this quantity a tiny value is the result of the difference between two terms that have completely different origins. Regarding the capabilities of detecting the axion dark matter, Fig. 3 shows the sensitivities of several experiments in the ma − gaγγ plane. In this plot, the thick blue and red lines are the regions where the axion abundance is responsible for all the observed DM. These lines were obtained by using jgj of order one. Moreover, the blue region corresponding to masses of the order of meV and gaγγ ≈ 10−12 GeV−1, lies very close to the projected IAXO sensitivity, so that it will be reachable in the near future. Looking back to our results we can conclude that this version of the 3 − 3 − 1 model, concerning the axion DM issue and the strong CP problem, is phenomenologically consistent. This model, besides its good qualities presented in the introduction, also possesses new degrees of freedom that are not yet experimentally probed. For instance, the model has charged and neutral scalars (besides the Higgs), extra vector bosons and extra quarks, that are expected to be heavy, and could, in principle, be searched at colliders. See Refs. [73,74] for recent studies concerning the 3 − 3 − 1 model phenomenology, in general, at the LHC. ACKNOWLEDGMENTS B. L. S. V. is thankful for the support of FAPESP funding Grant No. 2014/19164-6. A. R. R. C. would like to thank Coordenação de Aperfeiçoamento de Pessoal de Nível Superior (CAPES), Brasil, for financial support. J. C. MONTERO et al. PHYS. REV. D 97, 063015 (2018) 063015-10 [1] G. Bertone and D. Hooper, arXiv:1605.04909. [2] G. Bertone, D. Hooper, and J. Silk, Phys. Rep. 405, 279 (2005). [3] L. Bergström, Ann. Phys. (Berlin) 524, 479 (2012). [4] ATLAS-collaboration, J. High Energy Phys. 09 (2016) 175. [5] A. G. Dias, A. C. B. Machado, C. C. Nishi, A. Ringwald, and P. Vaudrevange, J. High Energy Phys. 06 (2014) 037. [6] P. W. Graham, I. G. Irastorza, S. K. Lamoreaux, A. Lindner, and K. A. van Bibber, Annu. Rev. Nucl. Part. Sci. 65, 485 (2015). [7] M. Battaglieri et al., arXiv:1707.04591. [8] J. E. Kim and G. Carosi, Rev. Mod. Phys. 82, 557 (2010). [9] J. Preskill, M. B. Wise, and F. Wilczek, Phys. Lett. 120B, 127 (1983). [10] L. Abbott and P. Sikivie, Phys. Lett. 120B, 133 (1983). [11] M. Dine and W. Fischler, Phys. Lett. 120B, 137 (1983). [12] C. D. R. Carvajal, B. L. Sánchez-Vega, and O. Zapata, Phys. Rev. D 96, 115035 (2017). [13] B. L. Sánchez-Vega and E. R. Schmitz, Phys. Rev. D 92, 053007 (2015). [14] B. L. Sánchez-Vega, J. C. Montero, and E. R. Schmitz, Phys. Rev. D 90, 055022 (2014). [15] C. A. d. S. Pires and O. P. Ravinez, Phys. Rev. D 58, 035008 (1998). [16] F. Pisano and V. Pleitez, Phys. Rev. D 46, 410 (1992). [17] P. H. Frampton, Phys. Rev. Lett. 69, 2889 (1992). [18] R. Foot, O. F. Hernández, F. Pisano, and V. Pleitez, Phys. Rev. D 47, 4158 (1993). [19] A. G. Dias, R. Martinez, and V. Pleitez, Eur. Phys. J. C 39, 101 (2005). [20] J. C. Montero and B. L. Sánchez-Vega, Phys. Rev. D 84, 055019 (2011). [21] R. Foot, H. N. Long, and T. A. Tran, Phys. Rev. D 50, R34 (1994). [22] J. C. Montero, F. Pisano, and V. Pleitez, Phys. Rev. D 47, 2918 (1993). [23] C. A. S. de Pires and P. S. Rodrigues da Silva, Eur. Phys. J. C 36, 397 (2004). [24] A. G. Dias, C. A. de S. Pires, and P. S. R. da Silva, Phys. Lett. B 628, 85 (2005). [25] P. V. Dong, T. T. Huong, D. T. Huong, and H. N. Long, Phys. Rev. D 74, 053003 (2006). [26] P. V. Dong, H. N. Long, and D. V. Soa, Phys. Rev. D 75, 073006 (2007). [27] J. K. Mizukoshi, C. A. de S. Pires, F. S. Queiroz, and P. S. R. da Silva, Phys. Rev. D 83, 065024 (2011). [28] A. G. Dias and V. Pleitez, Phys. Rev. D 69, 077702 (2004). [29] D. Cogollo, A. X. Gonzalez-Morales, F. S. Queiroz, and P. R. Teles, J. Cosmol. Astropart. Phys. 11 (2014) 002. [30] J. C. Montero and B. L. Sánchez-Vega, Phys. Rev. D 91, 037302 (2015). [31] S. M. Boucenna, J. W. F. Valle, and A. Vicente, Phys. Rev. D 92, 053001 (2015). [32] L. Clavelli and T. C. Yang, Phys. Rev. D 10, 658 (1974). [33] B. W. Lee and S. Weinberg, Phys. Rev. Lett. 38, 1237 (1977). [34] B. W. Lee and R. E. Shrock, Phys. Rev. D 17, 2410 (1978). [35] M. Singer, Phys. Rev. D 19, 296 (1979). [36] M. Singer, J. W. F. Valle, and J. Schechter, Phys. Rev. D 22, 738 (1980). [37] B. L. Sánchez-Vega, E. R. Schmitz, and J. C. Montero, Eur. Phys. J. C 78, 166 (2018). [38] J. E. Kim, Phys. Rev. Lett. 43, 103 (1979). [39] M. Dine, W. Fischler, and M. Srednicki, Phys. Lett. 104B, 199 (1981). [40] H. Georgi and M. B. Wise, Phys. Lett. 116B, 123 (1982). [41] M. Srednicki, Nucl. Phys. B260, 689 (1985). [42] W. A. Bardeen, R. Peccei, and T. Yanagida, Nucl. Phys. B279, 401 (1987). [43] M. Kamionkowski and J. March-Russell, Phys. Lett. B 282, 137 (1992). [44] R. Holman, S. D. Hsu, T. W. Kephart, E. W. Kolb, R. Watkins, and L. M. Widrow, Phys. Lett. B 282, 132 (1992). [45] L. M. Krauss and F. Wilczek, Phys. Rev. Lett. 62, 1221 (1989). [46] L. E. Ibáñez, Nucl. Phys. B398, 301 (1993). [47] K. S. Babu, I. Gogoladze, and K. Wang, Phys. Lett. B 560, 214 (2003). [48] K. Babu, I. Gogoladze, and K.Wang, Nucl. Phys.B660, 322 (2003). [49] L. E. Ibáñez and G. G. Ross, Phys. Lett. B 260, 291 (1991). [50] T. Banks and M. Dine, Phys. Rev. D 45, 1424 (1992). [51] C. Luhn and P. Ramond, J. High Energy Phys. 07 (2008) 085. [52] M. S. Turner and F. Wilczek, Phys. Rev. Lett. 66, 5 (1991). [53] D. Harari and P. Sikivie, Phys. Lett. B 195, 361 (1987). [54] A. Vilenkin and A. E. Everett, Phys. Rev. Lett. 48, 1867 (1982). [55] M. Kawasaki, K. Saikawa, and T. Sekiguchi, Phys. Rev. D 91, 065014 (2015). [56] O. Wantz and E. P. S. Shellard, Nucl. Phys. B829, 110 (2010). [57] O. Wantz and E. P. S. Shellard, Phys. Rev. D 82, 123508 (2010). [58] E. Kolb and M. Turner, The Early Universe (Addison- Wesley, New York, 1994), Vol. 69. [59] A. Vilenkin, Phys. Rep. 121, 263 (1985). [60] R. L. Davis, Phys. Lett. B 180, 225 (1986). [61] T. Hiramatsu, M. Kawasaki, T. Sekiguchi, M. Yamaguchi, and J. Yokoyama, Phys. Rev. D 83, 123531 (2011). [62] T. Hiramatsu, M. Kawasaki, K. Saikawa, and T. Sekiguchi, Phys. Rev. D 85, 105020 (2012). [63] P. Sikivie, Phys. Rev. Lett. 48, 1156 (1982). [64] A. Vilenkin and E. P. S. Shellard, Cosmic Strings and Other Topological Defects, Cambridge Monographs on Mathematical Physics (Cambridge University Press, Cambridge, England, 1994). [65] D. H. Lyth, Phys. Lett. B 275, 279 (1992). [66] A. Ringwald and K. Saikawa, Phys. Rev. D 93, 085031 (2016). [67] P. A. R. Ade et al. (Planck Collaboration), Astron. As- trophys. 594, A13 (2016). [68] C. D. R. Carvajal, A. G. Dias, C. C. Nishi, and B. L. Sánchez-Vega, J. High Energy Phys. 05 (2015) 069; 08 (2015) 103(E). [69] T. Dafni et al., Nuclear and Particle Physics Proceedings 273, 244 (2016). AXION DARK MATTER IN A 3 − 3 − 1 MODEL PHYS. REV. D 97, 063015 (2018) 063015-11 http://arXiv.org/abs/1605.04909 https://doi.org/10.1016/j.physrep.2004.08.031 https://doi.org/10.1016/j.physrep.2004.08.031 https://doi.org/10.1002/andp.201200116 https://doi.org/10.1007/JHEP09(2016)175 https://doi.org/10.1007/JHEP06(2014)037 https://doi.org/10.1146/annurev-nucl-102014-022120 https://doi.org/10.1146/annurev-nucl-102014-022120 http://arXiv.org/abs/1707.04591 https://doi.org/10.1103/RevModPhys.82.557 https://doi.org/10.1016/0370-2693(83)90637-8 https://doi.org/10.1016/0370-2693(83)90637-8 https://doi.org/10.1016/0370-2693(83)90638-X https://doi.org/10.1016/0370-2693(83)90639-1 https://doi.org/10.1103/PhysRevD.96.115035 https://doi.org/10.1103/PhysRevD.96.115035 https://doi.org/10.1103/PhysRevD.92.053007 https://doi.org/10.1103/PhysRevD.92.053007 https://doi.org/10.1103/PhysRevD.90.055022 https://doi.org/10.1103/PhysRevD.58.035008 https://doi.org/10.1103/PhysRevD.58.035008 https://doi.org/10.1103/PhysRevD.46.410 https://doi.org/10.1103/PhysRevLett.69.2889 https://doi.org/10.1103/PhysRevD.47.4158 https://doi.org/10.1103/PhysRevD.47.4158 https://doi.org/10.1140/epjc/s2004-02083-0 https://doi.org/10.1140/epjc/s2004-02083-0 https://doi.org/10.1103/PhysRevD.84.055019 https://doi.org/10.1103/PhysRevD.84.055019 https://doi.org/10.1103/PhysRevD.50.R34 https://doi.org/10.1103/PhysRevD.50.R34 https://doi.org/10.1103/PhysRevD.47.2918 https://doi.org/10.1103/PhysRevD.47.2918 https://doi.org/10.1140/epjc/s2004-01949-3 https://doi.org/10.1140/epjc/s2004-01949-3 https://doi.org/10.1016/j.physletb.2005.09.028 https://doi.org/10.1016/j.physletb.2005.09.028 https://doi.org/10.1103/PhysRevD.74.053003 https://doi.org/10.1103/PhysRevD.75.073006 https://doi.org/10.1103/PhysRevD.75.073006 https://doi.org/10.1103/PhysRevD.83.065024 https://doi.org/10.1103/PhysRevD.69.077702 https://doi.org/10.1103/PhysRevD.69.077702 https://doi.org/10.1088/1475-7516/2014/11/002 https://doi.org/10.1103/PhysRevD.91.037302 https://doi.org/10.1103/PhysRevD.91.037302 https://doi.org/10.1103/PhysRevD.92.053001 https://doi.org/10.1103/PhysRevD.92.053001 https://doi.org/10.1103/PhysRevD.10.658 https://doi.org/10.1103/PhysRevLett.38.1237 https://doi.org/10.1103/PhysRevLett.38.1237 https://doi.org/10.1103/PhysRevD.17.2410 https://doi.org/10.1103/PhysRevD.17.2410 https://doi.org/10.1103/PhysRevD.19.296 https://doi.org/10.1103/PhysRevD.22.738 https://doi.org/10.1103/PhysRevD.22.738 https://doi.org/10.1140/epjc/s10052-018-5626-2 https://doi.org/10.1103/PhysRevLett.43.103 https://doi.org/10.1016/0370-2693(81)90590-6 https://doi.org/10.1016/0370-2693(81)90590-6 https://doi.org/10.1016/0370-2693(82)90989-3 https://doi.org/10.1016/0550-3213(85)90054-9 https://doi.org/10.1016/0550-3213(87)90003-4 https://doi.org/10.1016/0550-3213(87)90003-4 https://doi.org/10.1016/0370-2693(92)90492-M https://doi.org/10.1016/0370-2693(92)90492-M https://doi.org/10.1016/0370-2693(92)90491-L https://doi.org/10.1103/PhysRevLett.62.1221 https://doi.org/10.1103/PhysRevLett.62.1221 https://doi.org/10.1016/0550-3213(93)90111-2 https://doi.org/10.1016/S0370-2693(03)00411-8 https://doi.org/10.1016/S0370-2693(03)00411-8 https://doi.org/10.1016/S0550-3213(03)00258-X https://doi.org/10.1016/S0550-3213(03)00258-X https://doi.org/10.1016/0370-2693(91)91614-2 https://doi.org/10.1016/0370-2693(91)91614-2 https://doi.org/10.1103/PhysRevD.45.1424 https://doi.org/10.1088/1126-6708/2008/07/085 https://doi.org/10.1103/PhysRevLett.66.5 https://doi.org/10.1103/PhysRevLett.66.5 https://doi.org/10.1016/0370-2693(87)90032-3 https://doi.org/10.1103/PhysRevLett.48.1867 https://doi.org/10.1103/PhysRevLett.48.1867 https://doi.org/10.1103/PhysRevD.91.065014 https://doi.org/10.1103/PhysRevD.91.065014 https://doi.org/10.1016/j.nuclphysb.2009.12.005 https://doi.org/10.1016/j.nuclphysb.2009.12.005 https://doi.org/10.1103/PhysRevD.82.123508 https://doi.org/10.1103/PhysRevD.82.123508 https://doi.org/10.1016/0370-1573(85)90033-X https://doi.org/10.1016/0370-2693(86)90300-X https://doi.org/10.1103/PhysRevD.83.123531 https://doi.org/10.1103/PhysRevD.85.105020 https://doi.org/10.1103/PhysRevLett.48.1156 https://doi.org/10.1016/0370-2693(92)91590-6 https://doi.org/10.1103/PhysRevD.93.085031 https://doi.org/10.1103/PhysRevD.93.085031 https://doi.org/10.1051/0004-6361/201525830 https://doi.org/10.1051/0004-6361/201525830 https://doi.org/10.1007/JHEP05(2015)069 https://doi.org/10.1007/JHEP08(2015)103 https://doi.org/10.1007/JHEP08(2015)103 https://doi.org/10.1016/j.nuclphysbps.2015.09.033 https://doi.org/10.1016/j.nuclphysbps.2015.09.033 [70] N. Bastidon (ALPS II Collaboration), in Proceedings, 11th PatrasWorkshoponAxions,WIMPsandWISPs (Axion-WIMP 2015): Zaragoza, Spain, June 22-26, 2015 (Verlag Deutsches Elektronen-Synchrotron, Hamburg, 2015), p. 31. [71] R. Bähre, B. Döbrich, J. Dreyling-Eschweiler, S. Ghazaryan, R. Hodajerdi, D. Horns, F. Januschek, E.-A. Knabbe, A. Lindner, D. Notz et al., J. Instrum. 8, T09001 (2013). [72] E. Armengaud, F. Avignone, M. Betz, P. Brax, P. Brun, G. Cantatore, J. Carmona, G. Carosi, F. Caspers, S. Caspi et al., J. Instrum. 9, T05002 (2014). [73] Q.-H. Cao and D.-M. Zhang, arXiv:1611.09337. [74] G. Corcella, C. Corianò, A. Costantini, and P. H. Frampton, Phys. Lett. B 773, 544 (2017). J. C. MONTERO et al. PHYS. REV. D 97, 063015 (2018) 063015-12 https://doi.org/10.1088/1748-0221/8/09/T09001 https://doi.org/10.1088/1748-0221/9/05/T05002 http://arXiv.org/abs/1611.09337 https://doi.org/10.1016/j.physletb.2017.09.015