New axion and hidden photon constraints from a solar data global fit

We present a new statistical analysis that combines helioseismology (sound speed, surface helium and convective radius) and solar neutrino observations (the $^8$B and $^7$Be fluxes) to place upper limits to the properties of non standard weakly interacting particles. Our analysis includes theoretical and observational errors, accounts for tensions between input parameters of solar models and can be easily extended to include other observational constraints. We present two applications to test the method: the well studied case of axions and axion-like particles and the more novel case of low mass hidden photons. For axions we obtain an upper limit at $3\sigma$ for the axion-photon coupling constant of $g_{a\gamma}\,<\,4.1 \cdot 10^{-10} \rm{GeV^{-1}}$. For hidden photons we obtain the most restrictive upper limit available accross a wide range of masses for the product of the kinetic mixing and mass of $\chi m<1.8 \cdot 10^{-12} \rm{eV}$ at $3\sigma$. Both cases improve the previous solar constraints based on the Standard Solar Models showing the power of using a global statistical approach.


Introduction and motivation
The Standard Model of Particles (SM) is not a complete theory, and extensions are necessary to address some of the most pressing open questions in fundamental physics. Among others, some of these are: nature of dark matter, matter-antimatter asymmetry in the Universe, origin of neutrino masses, strong CP-violation problem. In order to solve these problems, physics beyond the SM is needed and the existence of new particles and/or non-standard properties of known particles are generally invoked. Stellar physics, while far from being a closed subject, is a mature discipline and it provides us with an accurate understanding of the internal structure of stars and their evolution. The sheer mass and size, the extreme conditions of stellar matter and, in many cases, the extremely long lifetimes that allow to integrate small effects over a very long time make stars appealing as laboratories for particle physics under conditions not reproducible anywhere else.
In particular, the existence of weakly interacting particles beyond the standard model or of standard particles with non-standard properties can modify the internal structure and its evolution in different ways. For example, weakly interacting massive particles (WIMPs), too massive to be created inside stars, can be accreted from the dark matter halo of the Milky Way and contribute to the transport of energy inside stars [1][2][3] and, in the case of self-annihilating particles act as a localized energy source. Weakly interacting light particles, on the other hand, can be thermally produced in the stellar interiors and can easily escape due to their large mean free path and act as energy sinks. The increased rate of energy-loss then changes the internal structure of stars and also modifies evolutionary timescales. The production rate of each particle depends on the stellar conditions and, consequently, each will have an impact on different stars or on different evolutionary phases. Several works have given constraints on the properties of the weakly interacting light particles using different stellar objects. Some examples are: Globular clusters, where the luminosity of stars at the tip of the red giant branch [4][5][6] and the lifetime of horizontal branch stars [7][8][9][10] can be used to give constraints on axion (e.g. strength of coupling to photons and electrons) and neutrino (magnetic dipole moment) properties; massive stars, which can loose the blue loop thus contradicting the observation of Cepheids [11]; white dwarf stars, where the anomalous cooling induced by particles like axions shortens the cooling timescale of white dwarfs and impacts both on their luminosity function (LF) [12][13][14][15] and pulsation properties [16][17][18]; neutron stars, where additional axion energy losses [19,20] can even lead to observable deviations of surface cooling in real time [21]; SN1987A, where the extra cooling of the proto-neutron core would shorten the neutrino pulse [22][23][24][25] and the axions produced could produce a prompt gamma ray signal [23,26]. These ideas have been widely discussed in the literature being the most comprehensive reference the seminal book by G. Raffelt [27], which can be complemented with recent topical reviews on QCD axions [28,29], axion-like particles [30] and hidden photons [9].
Many studies have focused on using the Sun for setting limits on the properties of different types of particles; we review some of them below. The Sun is by far the bestknown star. The solar structure, revealed by helioseismology and solar neutrinos, is well determined, and accurate solar models give us information about the past, present and the future of he Sun [31]. While in some cases (e.g. axions) the most restrictive limits are not inferred from solar studies, the Sun remains the most useful benchmark for testing and validating both stellar models and, as it is partly the case of the present paper, different statistical approaches to constrain particle properties. Also, it is important to keep in mind that CAST [32] and the forthcoming IAXO [33,34] are experiments specifically designed to detect exotic particles directly from the Sun, so having predictions of expected solar fluxes for exotic particles remains an important aspect to consider. Solar constraints on particle properties have been generally derived from applying limits to variations of either neutrino fluxes [35,36] or the sound speed profile derived from helioseismology [35]. However, a systematic approach aimed at combining different sources of data accounting in detail for observation and theoretical errors is badly missing in the literature. Here, we try to supply such a tool.
The initial goal of this work is to extend the general statistical approach presented in [37] to constrain properties of particles (e.g. mass, coupling constant) making the best possible use of the available information on the Sun, both observational and theoretical. To do so we use the helioseismic data combined with the neutrino fluxes in a statistical approach that includes the theoretical and observational uncertainties and takes into account possible tensions among data and solar model input parameters. We then derive solar limits for the well-studied hadronic axions -to gauge the performance of our statistical approach-and for the more novel case of hidden photons for which the Sun sets the most restrictive limits on the kinetic mixing parameter for small hidden photon masses, m eV.
Axions are light pseudoscalar particles that were introduced by the Peccei-Quinn [38] solution to the strong CP problem. They are very light particles and interact with ordinary particles much like neutral mesons (π 0 , η,...) but with coupling strengths vastly weaker, and some model dependencies. The most relevant for astrophysics is whether axions couple with electrons with similar strength than to nucleons, because this coupling tends to be very efficient to produce axions in stellar environments. The DFSZ model [39] as any axion model embedded in a GUT are examples of axion models with tree-level axion coupling. Here, the main axion production mechanisms are the ABC proceses: axio-recombination, bremsstrahlung and Compton [40]. In the KSVZ model (hadronic axions) [41,42] this coupling is absent at tree level and the relevant axion coupling is to two photons. Here, the Primakoff effect is the mechanism that converts photons into an axions in the presence of electric or magnetic fields. The inverse Primakoff effect is used in helioscopes like CAST and SUMICO (and the proposed IAXO) to convert solar axions into detectable X-rays and constraint the solar axion flux.
The phenomenology of axions is extended to axion-like particles (ALPs) in a completely straight-forward way. ALPs are also bosons with a two-photon coupling, but their mass and interaction strength are in principle unrelated, unlike the axion case. These particles arise copiously in string compactifications (a string axiverse [43][44][45]) , where O(100) ALPs can be generically invoked, although the number of light enough particles to be relevant for stellar physics is completely model-dependent. ALPs can be behind some very interesting hints recently pointed out in astrophysics like the transparency of the universe to gamma-rays [46], the 3.55 keV line [47][48][49]and the soft excess of the coma cluster. [50] In this paper, we are mostly interested in the case of hadronic axions because the Sun is one of the most sensitive environments to look for their effects. Solar constraints are even more important for the case of ALPs, where the SN1987A constraint is absent in generic models (ALPs don't have generically large couplings to nucleons but axions do). We shall not endeavour constraining the axion-electron coupling because limits from white dwarves and red-giant stars are comparatively much stronger.
Until now, several works have given constraints to the axion-photon coupling constant g aγ using the variations that axions produce on helioseismologic quantities or neutrino fluxes. An upper limit for g aγ = 10 · 10 −10 GeV −1 is found by limiting the deviations of the solar sound speed profiles of solar models including axions [35]. The work of [35] also gives values for the solar neutrino fluxes depending on axion emission and [36] uses this relation to give a more restrictive constraint of g aγ = 7·10 −10 GeV −1 at a 3σ confidence level using the observed Φ( 8 B) solar neutrino flux measured by the SNO experiment (see [51] for the global analysis of the three SNO phases). In [52] they construct the so-called seismic models (non-standard solar models constructed ad-hoc to reproduce the sound speed derived from helioseismology), with different values for the axion-coupling constant and obtain an upper limit of g aγ = 2.5 · 10 −10 GeV −1 by comparing the predicted Φ( 8 B) with the experimental result with 1σ errors. For the mass range m a 0.02 eV, the most restrictive limit comes from the helioscope CAST with g aγ < 0.88 · 10 −10 GeV −1 [53]. The future helioscope IAXO should improve these results, as it is expected to reach sensitivities to the axion-coupling constant 1-1.5 orders of magnitude better than CAST [33].
Hidden photons (HPs), our second case of study, are vector bosons that couple weakly via kinetic mixing with standard photons. The kinetic mixing is represented by the parameter χ and together with the hidden photon mass, m, are the parameters that must be constrained (see [9,[54][55][56] and references therein). In fact, as described in Section 2.3, solar constraints are mostly sensitive to the product χm if the mass is below ∼ eV.
Hidden photons can only be produced from photon↔ oscillations, which are affected by the photon refraction in the solar plasma. The oscillations are resonant when the HP and photon dispersion relation match and this happens differently for transversely polarised photons and longitudinal excitations (plasmons). Resonant emission in the Sun is possible for HP masses below ∼0.3 keV, the highest plasma frequency ω p , which is also the highest photon effective mass. The emission of L-modes is more important for low HP masses (below m ∼eV) [8,9], for which resonant conditions happen all through the solar interior (each region of the Sun emits L-HPs with energy equal to the local plasma frequency). This is the case we will be interested in in this paper. Resonant emission of T-modes dominates the energy loss in the HP mass range ∼ eV-0.3 keV, and localises in a narrow spherical shell of the solar interior for which the HP mass matches the plasma frequency m ≃ ω p . This case is very interesting too, but our solar model codes need to be tuned for such highly localised sources, a task that we plan to endeavour in a forthcoming publication. For even higher masses, arguments from horizontal branch stars can give better constraints because, owing to the higher plasma frequencies, resonant production can happen up to higher masses, and when it is not possible the higher temperatures ensure the existence of photons of sufficient energy to account at least for the HP mass. Therefore, it is not crucial to reexamine this regime.
This paper is structured in the following way. In Section 2 we present the theoretical and experimental aspects of this work: solar models and data used and their respective treatment of uncertainties. In Section 3 we introduce the statistical method, in Section 4 we describe our main results and present new upper bounds to the axion mass and kinetic coupling of hidden photons. In Section 5 we compare our work to previous results and draw some conclusions and a short summary in Section 6.

Standard Solar Models
In this work we use standard solar models (SSM) as reference models. SSMs are calibrated to match the present-day solar radius R ⊙ = 6.9598·10 10 cm, luminosity L ⊙ = 3.8418·10 33 erg s −1 and surface metal-to-hydrogen ratio (Z/X) ⊙ . The choice of this last constraint is critical because it almost entirely determines the composition of the solar model and it has been the subject of much discussion over recent years in the context of the solar abundance problem [57][58][59][60].
One of the most used solar abundance composition is the GS98 compilation [61], based on 1-D solar atmosphere models, from which (Z/X) ⊙ = 0.0229. SSMs based on GS98 composition lead to very good agreement with the results from helioseismology and neutrino observables. The most recent solar abundance compilation by Asplund and collaborators [62] (hereafter AGSS09) results in much lower metallicity values with (Z/X) ⊙ = 0.0178. It uses a more accurate analysis with new three-dimensional hydrodynamical models of the solar atmosphere that reproduces very well the observed solar atmosphere (granulation, convective velocities and limb darkening among others) but it leads to SSMs with structures in strong disagreement with the solar structure inferred from helioseismology. SSM computed with GS98 presents a good agreement with the observed value of R CZ and Y S with 0.2σ and 1.1σ discrepancy while when using AGSS09 composition to build the SSM, the discrepancies are 2.5σ and 3.4σ respectively. For the sound speed profile the situation is similar with average rms values a factor 4 worse for the AGSS09 models compared with the GS98 ones [31], [37].
In [63], and more recently in [37], it was pointed out that current seismic observables are only sensitive to the actual opacity profile 1 . It is not possible with current helioseismic and solar neutrino data to constrain the solar composition independently of the radiative opacities used in solar models. Therefore, we cannot determine if the GS98 composition is correct, if there is a missing opacity source in current opacity calculations or if there is an additional mechanism of energy transport not accounted for in current solar models. However, for our purposes in this paper this is of little consequence.
The crucial aspect for our work here is that SSMs with GS98 composition reproduce very well the thermal stratification of the Sun. The production and effects of the exotic particles considered in this paper depend on the thermal structure of the Sun because energy-loss rates depend mainly on temperature and density, not on the detailed composition (see next section). Therefore, and regardless of the actual composition of the Sun, the SSM computed with the GS98 abundances offers a good reference model against which to test the impact of exotic particles on solar structure and properties. We have nevertheless considered models with both solar compositions for deriving more robust bounds on the properties of both axions and hidden photons. Models based on AGSS09 ones will let us also consider how these bounds depend on the reference model used.
We have computed all our models using GARSTEC (GARching STEllar Code) [65]. The physics used in the calculation of SSMs is as follows: OPAL05 as equation of state [66], OP opacities [67] complemented at low temperatures with molecular opacities [68] and nuclear rates from Solar Fusion II [69]. All models include microscopic diffusion as described in [70]. For more details about the SSMs used on this work see [31].

Models with axions
In the case of the Sun, production of axions occurs via the Primakoff effect, where axions are produced by the conversion of photons in the presence of the electric field of nuclei and electrons with an interaction lagrangian L aγ = g aγ B · Ea, where a is the axion field. Thus, constraints can be placed on axion-photon coupling constant g aγ . The corresponding energy-loss rate per unit mass [35] is given by equations 2.1-2.3.
where T is the temperature, ρ the density and F (κ 2 ) is a dimensionless function describing the screening effects. Here α is the fine structure constant, n B the baryon density, Y e the electrons per baryon, Y j = X j /A j and X j , A j and Z j represent, respectively, the mass fraction, atomic weight and atomic number of the nuclear species j. For solar conditions, [35] showed that the function F (κ 2 ) can be approximated by, We have computed two different solar models varying g aγ , each one using a different solar composition (GS98 or AGSS09). Each set is composed by 21 models with g aγ = g 10 · 10 −10 GeV −1 , where g 10 ranges from 0 to 20 with an interval of ∆g 10 = 1.

Models with hidden photons
The production of HPs in the interior of the stars can be seen as γ − γ ′ oscillations, described by the following lagrangian where A µν and B µν are the field strengths of the photon and HP field, A µ and B respectively. The energy loss rate per unit mass is based on the approximation presented in [9]. We only take into account the resonant emission of longitudinal HPs, which dominate (details in [9,54]), and gives where χ is the kinetic mixing parameter, m the mass of the hidden photon and ω P is the characteristic plasma frequency. Typical values in the solar center are ω P ∼ 0.3 keV and T ∼ 1 keV. By expanding the exponential then, it can be seen that, to first order, the temperature dependence of the energy-loss rate can be considered linear. Indeed, we should include a threshold factor 1 − (m/ω P ) 2 , but in the range we are interested, m eV, it is completely irrelevant. Note that now the rate depends only on the product χm.
We have thus computed models varying the product χm of the kinetic mixing parameter and the hidden photon mass. For hidden photons we have also computed two sets of models, for the GS98 and AGSS09 compositions, and each comprise nine models varying χm from 0 to 8 · 10 −12 eV with interval of ∆χm = 1 · 10 −12 eV.

Observables
The observables used for the study are the neutrino fluxes Φ( 8 B) and Φ( 7 Be), the surface helium abundance (Y s ), the convective radius (R CZ ) and the solar sound speed (c s ). For more clarity, we follow [37] and separate them into scalars (Φ( 8 B), Φ( 7 Be), Y s and R CZ ) and the solar sound speed for which we use its radial profile.
For scalars, the observational values and errors, together with the corresponding SSM results, are specified in Table 1. We adopt the same model errors as in [37].   [37]. The third and fourth columns summarize the observational and experimental values and errors. Neutrino fluxes are in 10 9 cm −2 s −1 for the Φ( 7 Be) and For the the radial profile of the sound speed we take 30 radial points r i in the range 0.05 < r i /R ⊙ < 0.80 and use the sound speed c i in those points as our observables. For each model considered here we have inverted the sound speed profile from helioseismic data. Inversion has been done using the SOLA inversion technique. Frequencies used are the BISON-13 dataset complemented with data of MDI, GOLF and IRIS. More details on both the frequency dataset and inversion technique are given in [74].

Uncertainties
There are two different error sources in this work: the fully correlated theoretical errors, induced by the uncertainties in the input parameters of solar models, and the observational errors associated to the helioseismic data and solar neutrino fluxes.
The theoretical errors are calculated by propagating the errors δI of the input parameters I. As in [37] these input parameters are the age of the Sun, diffusion coefficients, luminosity, opacity , astrophysical factors (S 11 , S 33 , S 34 , S 17 , S e7 , S 1,14 ) and in this work we have also used the surface abundances of the heavy elements relevant for the SSM construction (C, N, O, Ne, Mg, Si, S and Fe). For each of the observables Q we calculate the error contribution of each input parameter I. Finally, the total theoretical error for Q is found by combining all the input error contributions (eq. 2.6). (see [37] and [59] for details) where C Q,I represents the fractional variation of the observable Q when a fractional variation δI is applied to the input I. This coefficient is calculated as Values of the power-law exponents B Q,I and δI are given in [75].
The observational uncertainties are treated as uncorrelated and for the scalar quantities their uncertainties U Q correspond to the observational errors summarized in 1. For the solar sound speed, a relevant contribution to the uncertainty is provided by errors in the inversion of helioseismic data. Therefore, at each radial point r i where the first term is the experimental contribution coming from helioseimic frequencies, while the second and dominant term is the so-called "statistical" error in the inversion procedure estimated by [76]. As in all previous works, we assume there is no correlation between sound speed derived at different radial points.

Method and statistical procedure
The statistical approach is based upon constructing a χ 2 function with one degree of freedom given by the particle property we aim at constraining. We have built this function using 34 different observable inputs, the Φ( 8 B) and Φ( 7 Be) neutrino fluxes, Y s , R CZ and c i for 30 different value of r/R ⊙ where r/R ⊙ < 0.80. We exclude the sound speed profile for r/R ⊙ > 0.80 because the temperature gradient is adiabatic in that region and the sound speed profile is the same in all models and in agreement with helioseismic results. For values close to r/R ⊙ ∼ 1.0, the predicted sound speed deviates from the observations due to the surfaces effects not accounted for models. However, solar structure in this region is determined by the surface properties (mainly surface temperature and gravity) and thus not influenced by the properties of exotic particles.
The fractional difference between the theoretical and observational data δ Q is defined as where Q th and Q obs denote the model and observed values respectively for the observable Q. For each observable, δ Q is affected by uncorrelated errors from neutrino experiments and helioseismic data and by correlated errors from other sources such as opacities, age, nuclear astrophysical factors and abundance of elements.
Following the formulation of [77], the χ 2 function that accounts for the uncorrelated errors from the experimental data and the correlated ones induced by the model input parameters is given by The correction − I ξ I C Q,I is introduced to account for the effect of the correlated errors.
To normalize the effect of this correction a penalty I ξ 2 I is introduced to the χ 2 function. Each ξ I is a univariate Gaussian random variable, and it is a measure of the systematic error for each source I. The quadratic sum of all the ξ I that minimize χ 2 is the systematic pull of the system. Pulls give information about tensions between input parameters [77] in the models.
This method is equivalent to use the covariance matrix in order to compute the chisquared function. The main advantage of this method over the one using the covariance is that in this case is possible to study the individual contributions to the chi-squared function. For more details about the statistical method applied to solar models we refer the reader to [37].

Results
Equations 2.1 and 2.5 give the dependence of the energy loss rates for axions and hidden photons with temperature and density. Different dependences translate into different changes in the structure of the Sun. To facilitate discussion and interpretation of results we show in Figure 1 the energy-loss rate distribution as a function of the stellar radius for axions (red) and hidden photons (blue). The differences can be easily seen. Axions have an energy-loss rate distribution localized in the inner region of the Sun (r < 0.4R ⊙ ), whereas the energyloss rate distribution is much broader for hidden photons. The difference is mostly due to the temperature dependence of the energy-loss distribution 2.5 and 2.1. For axions this dependence is roughly proportional to T 6 , steepest than for the hidden photons for which, as previously discussed, the dependence is approximately linear under solar conditions. As a rule of thumb, we expect that observables will be affected in a different way for each case. For example, axions will have a major effect in the core regions of the sound speed profile while for models with the hidden photons the sound speed profile will be modified in all regions. Notice, however, that although energy losses can be localized in certain regions of the Sun, structural changes will be, to some degree, present all over the Sun. Local changes in certain quantities can derive in variations in the whole structure [78].  For axions, in Figure 2, we can see in the upper panels how the scalar observables Y s and R CZ change. Y s decreases with increasing values of the axion-coupling constant because it is almost perfectly correlated to the initial helium Y ini which decreases because it is necessary to have higher initial amounts of hydrogen so that the constraint imposed by matching L ⊙ can be fulfilled. The change in R CZ is quite small, as the plot shows, due to the to the fact that energy losses are localized towards the innermost regions of the Sun.

Scalar observables
The lower panels present results for the neutrino fluxes. These fluxes increase with increasing g 10 predominantly as a result of higher core temperatures and reach values wide out of the model and experimental errors. The Φ( 8 B) relative changes are more important than those of Φ( 7 Be) due to its stronger sensitivity to temperature. As a consequence, we expect that 8 B neutrino measurements give stronger constraint to g 10 . We note here that Φ( 8 B) alone has already been used to place constraints on g aγ [36,52]. It is important mentioning that for both the Φ( 8 B) and the Φ( 7 Be) fluxes, experimental results are excellent, with very small uncertainties. The constraining power of both fluxes, and particularly Φ( 8 B), is currently limited by uncertainties in solar models.
For the hidden photons case, Figure 3, the changes of the scalar quantities have the same qualitative behavior. The decrease in Y s and the increase in the neutrino fluxes are again a reflection of the adjustments necessary to satisfy the solar luminosity constraint: an increased initial hydrogen abundance and an increased core temperature. It is again not surprising that Φ( 8 B) shows the larger relative variations with respect to the SSM values. The case of R CZ is instructive for showing how the different types of particles induce different changes in the solar structure. Let us consider models for axions and hidden photons for which the Φ( 8 B) flux is similar, for example ∼ 15 · 10 6 cm 2 s −1 . In the case of axions, this corresponds to models that show ∆R CZ ≈ 0.003R ⊙ . For hidden photons, the equivalent model has ∆R CZ ≈ 0.006R ⊙ , a factor of two larger. Clearly, for a given change in the central  conditions, determined by the sensitive Φ( 8 B), hidden photons introduce larger changes that axions in the outer layers. This result is because they produce a broader distribution owing to their lower temperature dependence of the emission rate discussed previously.
We warn the reader, however, to refrain from directly comparing results of Figures 2 and 3 because they are shown as functions of parameters not related to each other in any way. We come back to comparing relative changes by axions or hidden photons in Section 5.

Solar sound speed profile
In Figure 4  Analyzing the effect of axions on the sound speed profile, we see that the larger effect occurs mainly in the inner region of the Sun. This is in line to the previous discussion on the behavior of R CZ , which showed little variation for axion models. This implies that the stronger constraints on g 10 are obtained from the measurements of the sound speed at r/R ⊙ < 0.35. On the other hand, for hidden photons the effect is noticeable in all the regions meaning that all the profile contribute to constrain χm.
Let us now compare solar models with different compositions. In the inner region, r/R ⊙ < 0.35, models with different composition are very close to each other. Thus, we expect the constraining power of solar models will not be affected by the choice of solar composition. In fact, in the next section we show differences being almost negligible. For the hidden photons models, the differences between compositions are more noticeable in the outer region because the SSMs corresponding to each of the compositions used already have a different behavior. However, even in this case results from models with different compositions show little variations, as we discuss in the next section.

Global results
In this section we present the results of our global analysis that combines all the observables using the method described in Section 3. The best fit values for g 10 and χm are determined by minimizing the χ 2 and the allowed region are defined by studying the ∆χ 2 = χ 2 − χ 2 m in distribution. Namely, the N σ allowed ranges for g 10 and χm are obtained by cutting at ∆χ 2 = N σ with N = 3. In general (but not always), the best fit model corresponds to SSMs,for which we obtain χ 2 /d.o.f = 72.5/33 and χ 2 /d.o.f = 47.3/33 for the AGSS09 and GS98 compositions respectively.
Axions -The best fit is obtained for g 10 = 0 (SSM) when using AGSS09 composition and for g 10 = 2 when using GS98 composition. The occurrence of a minimum of the χ 2 for a non-null value of g 10 can be related to the tension between the position of the R CZ and the bump of the solar speed profile right below the convective zone. In fact, it disappears when the more external solar speed data points are removed from the study, confirming this interpretation. It should be noted that the minimum at g 10 = 2 in the GS98 case is not statistically significant being the χ 2 substantially equivalent to those of the SSM.
In the left panel of Figure 5 we show the quantity N = ∆χ 2 (left axis) and the corresponding ∆χ 2 values (right axis) as a function of g 10 . We obtain g 10 < 3.93 for GS98 and a value of g 10 < 3.48 for AGSS09 by cutting at 3σ. We choose the more conservative result corresponding to the GS98 composition as our preferred bound for g 10 .
It is instructive trying to understand how much each individual constraint contributes to the global bound on g 10 . To this aim, we have calculated the χ 2 function separately for each observable. In the right panel of Figure 5 we represent these results for the GS98 composition. The values of N σ (left axis) and ∆χ 2 (right axis) are represented against g 10 . The most restrictive limit comes from the sound speed profile that, alone, sets the limit g 10 < 4.2 at 3σ. Φ( 8 B) also gives a strong bound, g 10 < 5.1. The Φ( 7 Be) neutrino flux gives g 10 < 6.4 while surface helium and convective boundary are the less restrictive observables. Notice that the individual χ 2 values are just indicative of the contribution to the final results. However, it should be kept in mind the global result is not merely the sum of all the individual contributions because observables are correlated.
Our results show that the sound speed profile is the most restrictive constraint. It is well known that a difference between the SSM sound speed profile and the solar one appears in the tachocline, the region just below the convective boundary, even for the models with GS98 composition (Figure 4). This is due to the inaccurate results of SSMs in modeling the composition gradient, gravitational settling, macroscopic motions or opacity among other possibilities in this region (e.g. [79][80][81]). Therefore, it is important to check that these differences do not affect our final result. We have repeated the global analysis excluding the region of r/R ⊙ > 0.6 so for this test 20 points of the solar speed profile are used. The exclusion of the external sound speed data points, as it expected, improves the quality of the fit being χ 2 SSM /d.o.f. = 16.2/23 for GS98 composition and χ 2 SSM /d.o.f. = 42.6/23 for AGSS09. However, because we derive the limits from ∆χ 2 , the bound we obtain for g 10 at 3σ does not significantly change, as shown in Figure 6. By excluding the tachocline region, the bounds are g 10 < 3.44 and g 10 < 3.6 for GS98 and AGSS09 models respectively, a modest different with respect to the case where the full sound speed profile is used. The fact that the axion effect is much important in the inner region of the sound-speed profile (see Figure 4) explains why the results depend only weakly on the inclusion of the bump.
Considering all previous results we assign a conservative 3σ upper limit of g aγ = 4 · 10 −10 GeV −1 for the axion-photon coupling constant. In Fig. 7 we plot our constrain in the context of other relevant astrophysical constraints for hadronic KSVZ axions and axion-like particles, together with the prospects of the future IAXO. Hidden photons -We show the results for models including hidden photons in Figure 8. The left panel shows, for both GS98 and AGSS09 sets of solar models, the ∆χ 2 function using all the observables and the corresponding N σ limits. The right panel illustrates the results for individual observables for the GS98 set of models. For both compositions, the best fit corresponds to the respective SSMs. The bounds at the 3σ CL are χm < 1.66 · 10 −12 eV and χm < 1.8 · 10 −12 eV for AGSS09 and GS98 respectively using the complete sound speed profile. If the region r > 0.6R ⊙ is excluded, then the limits are only marginally different.
The individual χ 2 contributions have similar behavior as in the axion case. The sound speed is the most restrictive observable, giving a constraint of χm < 2.15 · 10 −12 eV. The neutrinos fluxes also give strong constraints with values of χm < 2.62 · 10 −12 eV for the Φ( 8 B) flux and χm < 3.36 · 10 −12 eV for the Φ( 7 Be). The surface helium and the radius of the convective zone allow a wide range of parameter values.  Finally, we choose again our most conservative result to establish χm = 1.8 · 10 −12 eV as an upper bound to product of the kinetic mixing parameter and mass of hidden photons independently of the solar models employed. In Fig. 9 we see our new constraint in the mass-mixing plane together with other concurrent limits. Notably, the constraint is the most stringent in the range from 3 · 10 −5 eV to 8 eV. Figure 9: Constraints on hidden photons with mass m and kinetic mixing with photons χ. The global fit presented in this paper provides a constraint (Sun-L (global)) that improves over the solar constraint obtained only with the solar neutrino argument (Sun-L (ν)) and the direct detection constraint from XENON10 data [82].

Discussion
Here we have revised the solar limits on the axion-photon coupling constant and the product of the kinetic mixing and mass of the hidden photons by combining different solar constraints within a statistical approach that accounts for both experimental and theoretical errors.
Previous works on the subject have often relied on setting an upper limit to the total energy loss carried away from the Sun expressed as a fraction of the solar photon luminosity L ⊙ . It is then useful to present our results as a function of L i /L ⊙ , where L i = L a or L hp and denotes the luminosity carried away by axions or hidden photons respectively. Also, L i /L ⊙ sets a natural scale for comparing the effects of axion and hidden photon losses to each other.
In Figure 10, global results from the previous section are presented again (only for the GS98 sets of models) but now as a function of L i /L ⊙ , left panel shows N σ (or ∆χ 2 ) and the right panel shows the relative change of the neutrino fluxes with respect the SSMs values. We observe that the 3σ constraint for the hidden photon corresponds to a contribution on the solar luminosity of only about 2% while for the axions this value is around 3%. For comparison, in [36] the upper limit to axions is set by demanding L a < 0.1L ⊙ . We find that a much more restrictive limit can actually be imposed.
It is interesting to note that a smaller fraction of L ⊙ is required for the hidden photons than for axions to set a given confidence level. We can explain this by considering how the observables change in both cases. First, in Figure 4 we have seen that the changes on the solar sound speed profile are localized towards the center for the axion case while for the hidden photons the whole profile is affected. This introduces somewhat stronger constraints for hidden photons than for axions. Additionally, the neutrino fluxes show an analogous response, as we can see in the right panel of Figure 10 that illustrates the fractional change of the neutrino fluxes with respect to the SSMs for axions (red lines) and hidden photons (blue lines). For the same L i /L ⊙ value, changes in the neutrino fluxes are larger for the hidden photons models. It was shown by [36] that L a and Φ( 8 B) can be linked by an analytic relation of the form where α = 4.6, based on the older generation of SSMs computed by [35]. Our calculations yield α = 4.4 for axions, very close to the previous result. Interestingly, for hidden photons we find that the same functional form can be used but with a much steeper relation given by α = 5.7, as could be anticipated by results plotted in Figure 10. These results reinforce the importance of performing self-consistent solar model calculations to account for the effect of exotic particles. Also, they show that assuming a universality in the constraints that are imposed, e.g. a given L i /L ⊙ value, can lead to biased bounds on the properties of particles.
Solar limits on the coupling of hadronic axions to photons are not the most restrictive ones. However, this is a well-studied subject and it allows us to gauge the performance of our global statistical approach by comparing our findings to results from previous investigations.
The first work that placed a limit on the axion-photon coupling constant using helioseismology is [35], where a set of solar models was calculated including effects of axions self-consistently. The bound g 10 < 10 was found by restricting the axion luminosity at L a < 0.2L ⊙ . For this limit, the deviations the sound speed profile are 0.8% at R = 0.1R ⊙ .
More recently, [36] found g 10 < 7 by using the experimental Φ( 8 B) determination and, based upon its agreement with SSMs values, assuming that the actual solar flux cannot exceed the SSM prediction by more than 50% (roughly equivalent to a limit of 3σ model uncertainty). As mentioned before, this is comparable to setting the limit L a < 0.1L ⊙ .
A different approach has been followed in [52], where a seismic model is used to determine the effect of the axions. Seismic models are static models of the Sun which, by construction, reproduce the sound speed profile. They do not account for the evolutionary history which determines the chemical composition profile of the present Sun. Axions were included in the seismic models by adding the energy loss in the energy equation and modifying the central temperature and the abundances with respect the standard case so that the sound speed and solar luminosity are recovered 2 . The constraint on g 10 is derived from the experimental value of Φ( 8 B) and it is set to g 10 < 2.5 at a 1σ CL. To compare, from Figure  5 we find a limit at a 1σ CL around g 10 < 1.5.
In comparison to these works, our result g 10 < 4 to a 3σ CL is much more restricting. This is a direct consequence of our global approach of combining consistently helioseismic and solar neutrino constraints. As discussed in the previous section, the sound speed profile is the most restrictive observational constraint. If we restrict our analysis to using the Φ( 8 B) flux, we obtain g 10 < 5 and a corresponding L a < 0.05L ⊙ limit on the axion luminosity.
Hidden photons have a younger history than axions in the literature, and previous bounds on the kinetic mixing parameter based on solar models are limited to [8,9]. In the first case, a very conservative limit was derived by assuming that L hp < L ⊙ , leading to χm < 1.4 · 10 −11 eV. In the second, the more restrictive upper bound χm < 4 · 10 −12 eV for masses smaller than m < 0.3keV was derived from the condition L hp < 0.1L ⊙ and Equation 5.1. We have improved this limit by including the sound speed profile in the analysis and a more consistent treatment of uncertainties in our global approach. As discussed above, the limit χm < 1.8·10 −12 eV is obtained from models for which L hp < 0.02L ⊙ , a much smaller fraction than employed in previous works.
Using the limits on the parameters we have found for axions and hidden photons, we present upper limits for the respective fluxes on Earth expected in direct detection experiments such as CAST or IAXO. Using Eq. 15 from [53] and our limit g 10 = 4 we obtain Φ a ∼ 6.0 · 10 12 cm −1 · s −1 . During the CAST working years, the limiting flux have been found to be more restrictive for a wide range of axion masses. [53] find a limit for the axion-photon coupling constant of g 10 = 0.88 for m a 0.2eV corresponding to Φ a ∼ 2.9 · 10 11 cm −2 · s −1 and in [83] they find an upper limit of g 10 = 2.17 for the mass range of 0.02 < m a < 0.39eV corresponding to a solar flux of Φ a ∼ 1.8 · 10 12 cm −2 · s −1 . As already stated in the introduction, our limit does not improve previous ones but we can confirm that, at the CAST limits, it is not expected that axions would have an effect measurable with helioseismology and solar neutrinos.
On the other hand, for hidden photons we have found the most restrictive limit to date, improving the existing ones by about a factor of 2. Using Eqs 4.10 and 4.11 from [9], the upper limit for the flux on the Earth we derive is Φ HP ∼ 3.2 · 10 14 cm −2 · s −1 corresponding to χm < 1.8 · 10 −12 eV.
We remark that our bounds on axion-photon coupling and hidden photons kinetic mixing from are not affected from the ongoing solar abundance problem, a crisis brought about by the latest generation of spectroscopic determinations of the solar photosphere composition. Indeed, we systematically studied the role of the composition used in solar models for the constraints on g aγ and χm. Our limits are, to a good approximation, independent of it. We have nevertheless taken in all cases the most conservative limits. This also indicates that exotic energy losses, even assuming the broad radial distribution produced by hidden photons, cannot be advocated as a possible explanation of the disagreement between solar models with AGSS09 composition and helioseismic data.
There are different ways in which to move forward. Our global approach represents a qualitative step forward in combining different sources of solar data for using the Sun as a laboratory for particle physics. But further data is available which we have not considered. For example, frequencies of low-degree oscillations can be combined in an advantageous manner to gain further insight on the solar core (e.g. [84,85]). Using these data can enhance the constraining power of helioseismology. The caveat is that correlations among different helioseismic observables have to be taken into account. We have already acknowledged that even the radial profile of the sound speed presents correlations not yet been quantified in the literature. This is the next step we will pursue for improving our statistical approach.
It is important to notice that current uncertainties are dominated by model uncertainties. Moreover, there is not a unique dominant error source. For example, the 14% error in the theoretical Φ( 8 B) results from 8% error in S 17 , 7% from radiative opacities, 9% from solar composition. Thus, reducing modeling errors is a difficult task. A similar situation happens for the sound speed profile, where the solar composition and the radiative opacities play comparable roles.
On the other hand, the Φ(pp) and Φ(pep) fluxes have the interesting property that their theoretical errors are small, 0.6% and 1.2% respectively and therefore experimental determinations would place constraints on non-standard energy losses almost completely independent on solar modeling uncertainties. Recently [86], Borexino has for the first time directly measured the solar Φ(pp) with an observational error of ∼ 10% and the solar Φ(pep) [87] with an observational error of ∼ 20%. Due to these large experimental errors, we have not included these fluxes as observational inputs because they would not contribute to the derived limits. In a future generation of experiments, reaching a 1% error in such a measurement could provide a cross comparison between the solar photon and neutrino luminosities. Considering the case of hidden photons, for which the solar constraint is the most restrictive one over a wide range of photon masses, such a measurement would translate into a hard limit χm < 3 · 10 −12 eV, provided the measured value agrees with SSM results, and smaller if it does not. This limit, albeit larger than our estimate, has the attractive advantage that is practically independent of solar modeling.

Summary
We have extended the statistical approach presented [37] that combines in a consistent manner both helioseismic and solar neutrino data for using the Sun for particle physics studies.
As a test case we have applied the method to the well studied case of hadronic axions and showed that previous solar bounds on the axion-photon coupling can be improved with our methodology. We derive a strong upper limit of g aγ = 4 · 10 −10 GeV at a 3σ CL, almost a factor of two better than previous results. As a further application, we considered energy losses by hidden photons, an interesting case because, due to the low temperature dependence of the production rate, the Sun offers the most restrictive limits over a wide range of photon masses. For this case we have derived a new upper limit for their kinetic mixing parameter χm = 1.82 · 10 −12 eV for hidden photons in the mass range m 0.3 keV, more than a factor of 2 improvement over previous results [9] and better than the direct detection constraint of XENON10 [82].
By comparing in detail results from axions and hidden photons, we also conclude that relations between dark luminosities and observables such as Φ( 8 B), e.g. Equation 5.1, depend on the type of particle under consideration and should be employed consistently to avoid biasing results. Finally, we have checked that including hidden photons in solar models with AGSS09 composition degrades the agreement between models and helioseismology, and thus does not help in mitigating the solar abundance problem.
The approach is of course general and will be extended to account for correlations among observables, a fact hitherto neglected in the literature, that will allow us to include additional helioseismic constraints.