The Photon Content of the Neutron

In this work, we complete our CT18qed study with the neutron's photon parton distribution function (PDF), which is essential for the nucleus scattering phenomenology. Two methods, CT18lux and CT18qed, based on the LUXqed formalism and the DGLAP evolution, respectively, to determine the neutron's photon PDF have been presented. Various low-$Q^2$ non-perturbative variations have been carefully examined, which are treated as additional uncertainties on top of those induced by quark and gluon PDFs. The impacts of the momentum sum rule as well as isospin symmetry violation have been explored, and turn out to be negligible. A detailed comparison with other neutron's photon PDF sets has been performed, which shows a great improvement in the precision and a reasonable uncertainty estimation in our results. Finally, two phenomenological implications are demonstrated with photon-initiated processes: neutrino-nucleus $W$-boson production, which is important for the near-future TeV--PeV neutrino observations, and the axion-like particle production at a high-energy muon beam-dump experiment.


Introduction
Similar to the proton as a composite hadron, the neutron also contains its internal structure, mostly consisting of quarks and gluons as partons.The corresponding parton distribution functions (PDFs) can describe its internal structure up to very high precision, based on the modern strong interaction theory, Quantum Chromodynamics (QCD).A precise determination of the neutron's parton content plays a critical role in the hadron and particle physics frontier, which deepens our understanding of the fundamental parton dynamics [1].
Even as an electrically neutral particle, the neutron can still have photon content, due to its internal charged particles.Experimentally, the neutron's magnetic moment is measured to be µ n = −1.913µN [2], where µ N = eℏ/2m p is the nuclear magneton.In a more modern picture, the neutron's electromagnetic property can be described with the corresponding form factors, which are measured in various scattering experiments.See Ref. [3] for an overview of the most updated experiments and the global fitted results.The photon content induced by the electromagnetic form factors can be ascribed to the elastic component, as the neutron remains intact in a scattering process.The neutron's elastic photon can be determined through the so-called equivalent photon approximation (EPA) [4][5][6][7], which maps the electromagnetic form factors to the equivalent photon number or spectrum.This approach was adopted to include the elastic component of the proton's inclusive photon in the CT14qed PDF set [8], while the neutron's elastic photon was assumed to be zero regarding its zero electric charges.
Meanwhile, the neutron also has an inelastic photon component, which corresponds to inelastic scattering processes.In a deep-inelastic scattering (DIS), the internal quarks and gluons of the neutron can be resolved by a deep virtual photon.In such a picture, the photon becomes completely a parton of a neutron, which co-evolves with other partons (quarks and gluon) in terms of the DGLAP equations [9][10][11][12].The first quantum electrodynamics (QED) PDF set to include the photon as a parton of neutron came out as the MRST2004qed [13].The initial photon PDF there is parameterized with the radiation from the "primordial" up and down quarks, governed by the corresponding current or constituent quark masses.A small isospin violation at µ 0 was introduced to the valence quarks, with the corresponding amount fixed by momentum conservation.In such a way, the photon and other partons of neutron and proton were determined through a global QCD+QED analysis.With the isospin symmetric assumption at µ 0 , NNPDF fitted proton and neutron's photon content simultaneously through a global analysis of DIS and Drell-Yan data and released the corresponding PDF set as NNPDF2.3qed[14].Contemporaneously in the CT14qed PDFs [8], the proton's initial inelastic photon was constructed with the radiation from the valence quarks, with the initial momentum fraction ⟨xγ p ⟩(µ 2 0 ) determined with the ZEUS data on the isolated photon production [15].The neutron's PDFs are related to the proton ones with an isospin symmetry violation determined through ⟨xγ p ⟩(µ 2 0 ).The idea of EPA was extended to determine the inelastic photon as well in terms of the inelastic structure functions (SFs) [16][17][18].This approach was inherited by the LUXqed group [19,20], who established a rigorous framework based on the collinear factorization beyond the leading order by including proper MS matching terms.With the precisely measured elastic and inelastic structure functions, the proton's photon PDF was constrained up to a high precision at the level of 1-2%.Afterward, this LUXqed formalism [19,20] was adopted to determine the proton's initial photon PDF at µ 0 in the QCD+QED global analyses, including CT18qed [21], MMHT2015qed [22]/MSHT20qed [23], and NNPDF3.1luxQED1[24], which advanced our understanding of the proton's structure in the precision frontier [25].Based on the MRST2004qed proton-neutron isospin relation [13], both the MMHT2015qed [22] and MSHT20qed [23] analyses have released the neutron PDFs as well.
In this work, we extend our recent CT18qed study [21] to include the neutron PDFs as well, to provide a key input for many phenomenological applications related to nucleus scattering.Different from the CT14qed neutron PDFs [8], we will include the elastic photon component in terms of the electromagnetic form factors, similar to the treatment in the MMHT2015qed [22] and MSHT20qed [23] PDF sets.Following the CT18qed proton study [21], the neutron's photon PDF can be generated either through applying the LUXqed formalism at any scale or through evolving the LUXqed initialized photon up to higher scales, which we dub as CT18lux and CT18qed PDF sets, respectively.In the lower energy regime, the structure functions receive various nonpertubative contributions, such as higher-twist and target-mass corrections, which will result in different photon PDFs.We ascribed these variations as a part of the corresponding PDF uncertainties.
This paper is organized as follows.Sec. 2 lays out the neutron's photon PDF determination.More specifically, in Sec.2.1, we discuss the neutron's elastic photon component, followed by the CT18lux inelastic component in Sec.2.2, using the LUXqed formalism.In Sec.2.3, we develop the DGLAP methodology to consistently address the momentum sum rule and isospin violation in the mixed QCD and QED PDF evolution.A comprehensive comparison between the first and second generations of the neutron's photon PDFs is pre-sented in Sec.2.4.Two examples of the phenomenological implications are demonstrated in Sec. 3 and conclusions come afterward in Sec. 4.

The LUXqed photon PDF formalism
As implied by the equivalent photon approximation [7], the nucleon's photon PDF is related to the corresponding structure functions.The elastic structure functions can be constructed in terms of the electromagnetic form factors as [7] where m is the nucleon mass.By mapping the equivalent photon spectrum to the elastic2 photon PDF, we obtain [7] xγ el = α π where the photon virtuality Q 2 can be related to the transverse momentum Q ⊥ as An improved form based on the Q 2 ⊥ integration including the proper integral interval can be found in Refs.[27,28].
This idea was applied to the inelastic photon as well in Refs.[16][17][18].With a rigorous collinear factorization framework, the LUXqed group extended this approach beyond the leading order by including a proper MS matching term [19,20], which we dub as the LUXqed formalism (2.4) ]/z corresponds to the leading order DGLAP splitting kernel.The α ph is the physical fine structure constant, defined as where Π(q 2 , µ 2 ) is the vacuum polarization function, evaluated at q 2 = −Q 2 .The first term in Eq. (2.4) with the square bracket is referred to as the physical factorization component, and the second one that only involves F 2 is the MS conversion [20].To perform the integration in Eq. (2.4), the neutron's structure functions F n 2,L in the complete (x, Q 2 ) plane are needed.With a precise knowledge of the structure functions F n 2,L , the neutron's photon PDF can be determined up to a high precision, which is the main goal of this work.Following the LUXqed treatment [19,20], we divide the (x, Q 2 ) plane mainly into three regions, with a smooth transition when crossing a boundary, illustrated in Fig. 1.In the high- PDF , perturbative QCD (pQCD) is well-established, thus we employ pQCD to calculate the structure functions.Here Q 2 PDF denotes the matching scale, which varies between 5 and 9 GeV 2 in this work, with the difference quantifying the corresponding uncertainty.When Q 2 decreases, we gradually enter the non-perturbative regions.The region with In this region, we directly take the structure functions measured by the HERMES experiment with the GD-11P and GD-11D fits [29] with the ALLM parametrization [30,31].Considering the potential loss of reliability in the extrapolated low-x region, x ≲ 10 −4 , we take the pQCD structure functions as an alternative choice.In the resonance region W 2 < W 2 res , the nucleon structure functions show many resonance features.In this work, we take the CLAS [32] and Christy-Bosted [33][34][35] fits to determine the resonance structure functions.Below the resonance threshold W 2 < (m + m π ) 2 , the nucleon inelastic structure functions vanish, which does not contribute.Finally, the elastic form factors are taken from the global fit of the elastic scattering data [3] or the Galster parameterization [36], shown as the red band around x ∼ 1 as the nucleon remains intact.A detailed description of the structure-function determination will come in the next two subsections.

The elastic photon PDF
In this subsection, we explore the nucleon's elastic photon, with the elastic form factors from fitting the global data [3].We will also compare our results with the MMHT15qed [22] and MSHT20qed [23] ones, which uses the Galster parameterization [36].
The LUXqed formalism, Eq. (2.4), can be applied to the elastic photon as well.The elastic structure functions in Eq. (2.1) can be rewritten as [37] where τ = Q 2 /4m 2 .Here δ(1 − x) factor indicates that the nucleon remains intact, which corresponds to the red narrow band in Fig. 1.Substituting Eq. (2.6) into Eq.(2.4) ends up with the LUXqed elastic photon as Here this LUXqed form shares the same philosophy as the EPA in Eq. (2.2), but improves with a proper matching, which facilitates the higher-order extension.Here the physical factorization and MS terms are merged together, resulting in the Q 2 integration up to infinity.We remind the readers that this upper integration limit is different from the one in the MMHT15qed [22] and MSHT20qed [23] fits, which adopted the starting scale µ 2 0 , while the elastic photon at higher scale was obtained through a derived evolution equation.See Ref. [22] for the details.
As mentioned in Sec. 1, the electromagnetic form factors G E,M (Q 2 ) can be directly extracted from the experimental measurements.In comparison with the proton form factors, the neutron ones are poorly known, as no free and stable neutron exists in a natural environment.Therefore, the neutron form factors have to be extracted from nuclei (such as Deuterium or Helium), by subtracting the dominant proton contribution after properly accounting for nuclear effects.See Ref. [3] for the details of the determination methodology.
Instead, we adopt the global fitted G n E,M results by Ye et al. [3] to determine the neutron's elastic photon component in this work.The comparison between the two different elastic form factors is displayed in Fig. 2.Here the global fitted G n E,M are displayed as the red curves, and the red thin bands ("par.")indicate the corresponding uncertainties parameterized with analytical formulas.The green dashed curves ("cov.")refer to the same uncertainties obtained with the covariance matrix of the fit, which agree with the parameterized forms very well through the whole Q 2 range.In practice, we adopt the analytical form in our calculation for convenience.In comparison with the Ye et al.'s fit, the Galster parameterization [36,38] captures the general feature of the neutron's form factors, while with an underestimation of the uncertainties, especially at a high photon virtuality Q 2 .
In Fig. 3, we show the neutron's elastic photon PDF, xγ el (x, µ 2 ), at µ = 1.3 GeV.The higher-scale photon PDFs behave similarly, with only a small variation resulting from the running of α(µ 2 ) in Eq. (2.7), as shown in the lower right panel.First, we see that the neutron's elastic photon is dominated by the magnetic contribution, while the electric Figure 3. Upper: the neutron's elastic photon at scale µ = 1.3 GeV, based on the Galster [36] and Ye [3]'s electromagnetic form factors. Lower: the comparison between neutron and proton's elastic photon PDFs, generated with the Ye's form factors.We also include MSHT20qed result [23] for a comparison.
form factor only contributes a percent level to the total elastic photon.This can be easily understood in terms of its zero electric charge, which is different from the proton case, as shown in the lower left panel of Fig. 3.In the proton case, the electric (magnetic) contribution dominates the elastic photon at low (large) x, with equality around x ∼ 0.2.Second, we also see that the neutron's total elastic photon PDF is generally much smaller than the proton's, by about two orders of magnitude at a small x.At a large x (x ≳ 0.7) where the magnetic contribution dominates, the neutron and proton elastic photons scale as the magnetic moment ratio µ n /µ p .Third, the neutron's elastic photon induced by the Galster parameterization is slightly larger than the one with Ye's form factors by about 3% at a low x, mainly due to the larger low-scale G M when Q 2 ≲ 0.6 GeV 2 as shown in Fig. 2. As we observed before, the Galster parameterization underestimated the uncertainties, which resulted in smaller error bands.Meanwhile, we include the elastic photon from the MSHT20qed fit [23] for a comparison, which is about 30% of our result, mainly resulting from a much lower integration limit, as we mentioned above.However, this difference of the proton's elastic photon between the one with Ye's form factors and the MSHT20qed result is much smaller, which only shows up when x ≳ 0.7 at µ = 1.3 GeV3 as shown in Fig. 3 lower left panel.

The inelastic photon PDF with the LUXqed approach: CT18lux
As discussed before, when applying the LUXqed formalism to the neutron's inelastic photon, we need to integrate Eq. (2.4) throughout the complete (x, Q 2 ) plane as shown in Fig. 1.
Here we take similar breakups as the original LUXqed treatment [20].That is, we divide the (x, Q 2 ) plane into regions discussed as follows.
Resonance.W 2 < W 2 lo = 3 GeV 2 is defined as the resonance region.We remind the readers that the hadronic inelastic kinematics is bounded from below with a threshold W 2 ≥ (m + m π ) 2 , shown as the lower boundary in Fig. 1.Similar to the LUXqed treatment, we take the neutron's resonance structure functions either from the CLAS experiment [32] or the Christy-Bosted fit in 2007 (CB07) [33,34] as well as the updated one in 2021 (CB21) [35].With a few representative scales Q 2 , we compare F 1,2 from these three resources in Fig. 4. We see an overall agreement along the whole W 2 range, even though the difference is slightly larger than the proton case, which can be found in our recent work [21].We take the CLAS fit as our default choice, while the variation to CB21 fits quantifies the resonance uncertainty.The neutron structure functions in the resonance region, taken from the CLAS experiment [32] and the Christy-Bosted fit [33][34][35], respectively.The "CH" denotes the smooth transition from the CLAS to HERMES with Eq. (2.17).
Low-Q 2 continuum.W 2 > W 2 hi = 4 GeV 2 and Q 2 < Q 2 lo is the low-Q 2 continuum region.We follow the LUXqed treatment [20] and construct the structure functions based on the HERMES GD11-P and GD11-D fits [29], which adopted a wide range of data and the ALLM functional form [30,31].The HERMES Collaboration measured the inclusive DIS cross section from both proton and deuteron targets.In Fig. 5, we show the sum of longitudinal and transverse photon-absorption cross sections σ T +L = σ T + σ L , with which we can directly construct the F 2 structure function as (2.12) In order to get the F 1 or F L , we need the longitudinal-to-transverse cross-section ratio the LHAPDF interpolation of its x grid.
which can be obtained from the R1998 fit [39] or Ref. [37].In such a way, the longitudinal structure function can be constructed as . (2.14) Here we always keep the target-mass factor (1 + 4x 2 m 2 /Q 2 ), which will be explored in more detail for its corrections.In order to obtain the neutron's structure functions, we take the assumption that by ignoring the Deuterium's nuclear corrections, with the result shown in Fig. 5 lower panel.We remind that σ d T +L is defined as a cross section per nucleon.The uncertainties induced by this low-Q 2 HERMES SFs are propagated from cross section quantified by The difference of R L/T ratio between Deuterium and proton ∆R = R d − R p is measured by many experiments, such as SLAC [40,41], and NMC [42][43][44].It is found that ∆R is compatible with zero within its uncertainty with size O(≲ 10%), through a wide x range, such as 0.002 ∼ 0.4 [44] and even up to x < 0.8 [40] .In this work, we take , with the uncertainty captured by the variation of R1998 [39] by ±50%.
W 2 transition.Following the LUXqed treatment [19,20], the structure functions in the gap between the resonance and low-Q 2 continuum regions (W 2 lo < W 2 < W 2 hi ) are dealt with a smooth transition that where a = 2, L and (2.18) GeV 2 is dubbed as the high-Q 2 continuum region.In this region, the perturbative quantum chromodynamics (pQCD) is well established, which we can rely on to calculate the corresponding structure functions up to a high order [45].Based on the CT18lux methodology [21], we need to know the quark and gluon PDFs at the scale when µ > µ 0 .In this work, we adopt the isospin symmetric approximation (ISA) to relate the proton and neutron quark PDFs, (−) Due to the QED evolution effect, a small isospin violation will be generated through the γ → q q splitting [8, 13, 14], which will be properly addressed later with the CT18qed approach [21].Similar to the LUXqed treatment [20], a variation as Q 2 PDF = 5 GeV 2 is introduced to explore the corresponding PDF matching uncertainty.
Q 2 transition.Similar to the low-W 2 transition region, we also smoothly transit the SFs from the low-Q 2 HERMES to high-Q 2 pQCD continuum region when The smooth structure functions are constructed similarly as where The structure functions F 2 , F L in this transition region with a few representative Bjorken-x's are shown in Fig. 6.Different from the low-W 2 transition as the "CH" shows in Fig. 4, we see sizable mismatches and error bands from the low-Q 2 HERMES to the high-Q 2 pQCD continuum region, especially at a small x.It mainly comes from the low-x extrapolation of the HERMES fit, reflecting the lack of data in this region.
Here in Fig. 6, we also include the SF uncertainty at high Q 2 induced by the pQCD calculation.We see that in the moderate-x region, such as 10 −3 ≲ x ≲ 0.5, the HERMES results agree with the pQCD ones well in the matching region.However, with x decreasing into the extrapolation region, x ≲ 10 −4 , we see a clear deviation, as the HERMES GD11-D fit only parameterize the deuteron SFs F d 2 in the region 0.89 × 10 −3 < x < 0.9 [29].It is questionable which prediction is more reliable in this small-x and small-Q 2 region (denoted as "sxQ2" later) between the HERMES extrapolation and the pQCD calculation.In this work, we take the HERMES fits as default to incorporate low-Q 2 SF data.Meanwhile, we also take the pQCD SFs as an alternative choice in this region (x ≲ 10 −4 ), with the difference as an additional quantification of the corresponding uncertainty.
The error bands opening to the low-Q 2 direction are taken from the HERMES data with Eq. (2.16), while the high-Q 2 bands indicate the pQCD PDF uncertainty.
At this stage, the inelastic photon PDF can be directly calculated with the LUXqed formalism, i.e., the CT18lux approach.In Fig. 7, we show the inelastic as well as the total photon PDFs in comparison with the elastic components at a few representative scales, µ = 1.3, 3, and 100 GeV.Similar to the proton case, we see the inelastic photon grows with scale very fast, especially in the small-x region.In the right panel of Fig. 7, we display the fractions of elastic and inelastic components.We see that due to the neutral electric charge, the elastic photon only contributes a small fraction to the total photon PDF within x ≲ O(a few • 10 −1 ), which differs from the proton case where the elastic photon takes over at a low scale when µ ≲ 10 GeV. 4 (For this reason, the CT14qed took a zero elastic component for neutron's photon PDF [8].)Only at a very large x, the neutron's elastic photon can play an important role, mostly induced by the neutron's magnetic form factor.We also notice a numerical subtlety here.As shown in Fig. 7 (right), we see that the inelastic component of the neutron's photon at µ = 1.3 GeV becomes zero around x ∼ 0.6, and becomes positive again when x ≳ 0.8.It is resulted from the subtracted structure function F 2 (x/z, µ 2 ) in Eq. (2.4) falls into the non-perturbative resonance region, which exceeds the Q 2 integration.However, considering the numerical absolute smallness (≲ 10 −6 ) and phenomenological irrelevance of the photon PDF at such a large x value, we just enforce it as zero in this region.The specific x range depends on the scale µ as well as other low-Q 2 SF treatments as well.
In the high-Q 2 continuum region as shown in Fig. 1, the theoretical uncertainties for the pQCD structure functions propagate into the inelastic photon, including the higher-twist (HT) as well as target-mass (TM) corrections.
The HT correction is explained as follows.In terms of the spirit of operator product expansion [46], the DIS hadronic tensor can be expanded with coefficient functions together with local operator matrix elements, where ν = q • p/m = E − E ′ is the lepton's energy loss and q 2 = −Q 2 is the photon virtuality.
The fermionic operators can be written as while the scalar and vector operators can be defined similarly [47].The operator twist τ = d − j is defined as the difference between the mass dimension d and the corresponding spin j [48].In the large virtuality limit Q 2 → ∞, the forward amplitude T µν (q 2 , ν) can be approximated as where x = Q 2 /(2mν) is the Bjorken scaling.The Wilson coefficients scale as where γ j is the anomalous dimension of operator O (j) i , which increase monotonically with j.In perturbative QCD, the leading-twist operators correspond to τ = 2 with γ 2 = 0 which are defined as parton distribution functions, while higher twists starting from τ = 4 serve as power corrections O(Λ 2 /Q 2 ).
In this work, we take the same treatment as the CT18qed study [21].The HT corrections to the F 2 are taken from the CJ15 NLO fitting [49], where GeV 2 with parameters h 0,1,2 fitted from data.The higher twist contribution to the longitudinal SF is taken as where A HT = 5.5 ± 0.6 GeV 2 from the HERA fit [50].The corresponding HT corrections to SFs, defined as are shown in Fig. 8 (left).Also, we see the HT generally gives larger corrections to F L than F 2 , mainly resulting from the smaller magnitude of F L with respect to that of F 2 .Besides the higher-twist corrections, the non-zero target mass (TM) of the scattered nucleus also introduces kinematic power corrections O(x 2 m 2 /Q 2 ) [51], which was first calculated by Georgi and Politzer [52].With the Nachtmann variable [53] we have the TM corrected SFs as [54,55] where and 2,L are the SFs in the m 2 /Q 2 → 0 limit.The TM corrected SFs normalized to the uncorrected ones as are shown in Fig. 8 (right).Similar to the HT case, we see the TM corrections only become significant at large x and small Q 2 .In Fig. 9, we show various resources contributing to the uncertainties of neutron's photon PDF (summing over elastic and inelastic components), based on the CT18lux methodology.Most of them are in parallel with the proton case studied in Ref. [21].We discuss them in sequence as follows.
• Elastic.The elastic photon variation is deduced from the uncertainty of electromagnetic form factors obtained in the Ye et al.'s fit [3].Due to a small contribution to the total photon, we see the elastic uncertainty is small in most of the x range, which is only noticeable up to a percent level at a very large momentum fraction x ≳ 0.5.• CB21.The "CB21" denotes the variation from our default CLAS resonance structure functions [32] to the CB21 [35] fit as shown in Fig. 4. Different from the proton case, we see that the CB21 resonance dominates the neutron's photon uncertainty, suggesting the imprecise determination of the neutron structure functions in the resonance region.
We have made sure that the CB07 [33,34] fit gives a similar size of variation, as shown in Fig. 10.Indeed, the neutron's SFs are not directly measured from experiments, which are instead extracted from the scatterings of Deuterium/Helium and proton targets in this region, based on some assumptions of the nuclear corrections.See Ref. [34] for an example.For this reason, the neutron's SFs are not as well determined as the proton's, which introduces larger uncertainty on the inelastic photon PDF.In addition, we see that with an increment of the scale µ, the CB21/07 resonance uncertainty decreases, because of the increasing contribution from the pQCD continuum region.Figure 10.The variation of the CT18lux neutron's photon PDF at µ = 1.3, 3, 100 GeV with the resonance structure functions of CB07 [33,34], CB21 [35] and CLAS [32] fits.
• HERMES.The neutron's SFs in the low-Q 2 HERMES region are constructed in terms of Eq. (2.15).In this work, the nuclear corrections are neglected in the central relation, while the uncertainty is conservatively propagated with the squared sum according to Eq. (2.16).We see that the HERMES uncertainty can induce up to an 8% variation for the neutron's photon in the small-x and small-µ region.The dying out HERMES' impact in the large x direction is due to the better matching of SFs in the transition region as shown in Fig. 6.With increasing scale µ, the HERMES-induced uncertainty gets the most pronounced around µ = 3 GeV.When µ > 3 GeV, this impact decreases due to a smaller contribution from the HERMES region, similar to the CB21 resonance impact.• R1998.As described above, we have assigned a large variation to the longitudinal-totransverse cross-section ratio as R L/T = R 1998 (1 ± 50%), which introduces the corresponding inelastic photon PDF uncertainty.This variation is as small as about 1% ∼ 2% relevant in the large x region when x > 0.1, due to the small contribution from the longitudinal SF F L in the HERMES region.• sxQ2.As shown in Fig. 6, the low-Q 2 SFs of the HERMES fit deviates from the pQCD ones in the extrapolated small-x region, x ≲ 10 −4 .In this work, we take the HERMES fit as our default choice, while the pQCD SFs as an alternative choice gives an additional error set to quantify the SF variation in this small-x and small-Q 2 region (shortened as "sxQ2").In terms of Fig. 9, the "sxQ2" SF variation gives a more than 10% variation for the photon PDF at a low x (x ≲ 10 −4 ) and a low scale (µ = 1 ∼ 3 GeV), which is the largest error source in this region.However, with the increasing of the scale, e.g., at µ = 100 GeV, the uncertainty from the sxQ2 SF variation falls quickly into the error band from other resources, as a result of the increasing contribution from the high-Q 2 pQCD region.• HT.The higher-twist corrections are described above for F HT 2,L .The corresponding uncertainty is negligible in most scenarios, as a result of the 1/Q 2 suppression in the pQCD region.
• Q2PDF.We take Q 2 PDF (flatten as "Q2PDF") to denote the matching scale between the low-Q 2 continuum HERMES and high-Q 2 pQCD regions.The matching uncertainty is quantified by varying the default Q 2 PDF = 9 GeV 2 to 5 GeV 2 .We see in both µ = 1.3 GeV and 100 GeV cases, the matching uncertainty is under control, while it can introduce a large variation up to 4% at x ∼ 10 −5 when µ = 3 GeV.At µ = 1.3(100)GeV, the contribution to the photon PDF from the pQCD (low-Q 2 HERMES) SFs is small.As a result, the variation of Q 2 PDF = 9(5) GeV 2 does not make a big difference in both cases.In contrast at µ = 3 GeV, where the photon PDF has a large contribution from the matching region.For this reason, the mismatch between the pQCD and HERMES SFs will be reflected in the photon PDF variation.The larger photon PDF with Q 2 PDF = 5 GeV 2 in the small x region is due to the larger low-Q 2 SFs in the HERMES extrapolation as shown in Fig. 6.In comparison, the proton's photon uncertainty from Q 2 PDF is much smaller [21], resulted from a much better matching for the proton SFs in these two regions.• MHO.The missing higher order (MHO) uncertainty is quantified by varying the upper integration limit µ 2 /(1 − x) in Eq. (2.4) to be µ 2 , with the corresponding MS matching term [20,21].We see that the MHO uncertainty only becomes significant in the x → 1 limit, as a result of the asymptotic divergence of µ 2 /(1 − x).Nevertheless, the MHO only gives 1% ∼ 2% uncertainty in the normal x region.Compared with µ = 1.3 and 100 GeV cases, the MHO variation of γ(x, 3 GeV) is slightly larger, due to the amplification of the transition impact in the • TM.In the small-x region, the target-mass corrections give a negligible effect, which only significantly deviates when x → 1.In addition, we see the TM effect increases with scale µ, which can even dominate when x > 0.2 at µ = 100 GeV.A similar behavior occurs in the proton case already but at a much larger momentum fraction around x > 0.6 [21].• q, g PDF unc.In Fig. 9, the gray error bands denote the uncertainty induced by the quark/gluon PDF error sets.At µ = 1.3 GeV, we see the q, g induced uncertainty is negligible, suggesting a small contribution from pQCD SFs.In contrast, the q, g induced uncertainty increases with scale, as a result of the increment of pQCD significance.At µ = 100 GeV, the q, g induced uncertainty is about 1% ∼ 3%, the same level as the proton one [21].To summarize this subsection, based on the CT18lux approach we present the neutron's photon PDFs (including the elastic and inelastic components) in Fig. 7. Various photon PDF uncertainties are shown in Fig. 9, with the details discussed above.Many variations share similarities with the proton case, while the impacts of the low-Q 2 CLAS/CB21 resonance and HERMES get enhanced, due to the larger uncertainties of neutron SFs.The neutron's elastic uncertainty is much smaller than the proton's, thanks to its smaller contribution.The HERMES and matching scale Q 2 PDF uncertainties get pronounced around µ = 3 GeV, resulting from the mismatch between the low-Q 2 HERMES and high-Q 2 pQCD continuum regions.

The DGLAP evolution approach: CT18qed
Similar to our previous study, the neutron's photon PDF can be determined through the DGLAP evolution, with the initialization at µ 0 taken from the LUXqed photon, which we dub as the "CT18qed" methodology [21].Momentum sum rule.As emphasized in Ref. [21], the additional photon in the CT18lux approach introduces a violation of the momentum sum rule, as where Σ = i (q i + qi ) is the flavor singlet.At the initialization scale µ 0 = 1.3 GeV, the CT18lux neutron's photon momentum fraction for the elastic and inelastic components are with the total listed in Table 1.In comparison with the proton case, Eq. ( 28) of Ref. [21], the neutron elastic photon is significantly smaller due to its neutral electric charge.The inelastic component is comparable to the proton one, with a slightly smaller size due to its smaller charge-weighted singlet Σ e = i e 2 i (q i + qi ), where e i is the charge of quark q i .We take out a momentum fraction from the gluon PDF to enforce the momentum sum rule, In such a way, the initial gluon PDF is re-scaled by a factor of Note that there exists an ambiguity about whether the momentum sum rule should include the elastic photon.Different from our previous treatment in the proton's QED PDFs [21], we only include the inelastic component in Eq. (2.35) based on two arguments.First, in the elastic scattering, the nucleon remains intact and therefore the elastic photon does not become an internal parton.Moreover, the DGLAP equations conserve the momentum sum rule in Eq. (2.35), which only involves the inelastic photon in evolution.

Isospin symmetry violation (ISV).
In the DGLAP QCD+QED co-evolution, the quark and gluon PDFs will receive modification from the photon splittings with respect to the pure QCD evolution.In the CT18lux calculation of the neutron's photon PDF, we have to rely on the assumption of the isospin symmetry, Eq. (2. 19), to obtain the quark and gluon PDFs at all scales.However, the isospin symmetry will be violated starting from the O(α) order, due to different electric charges of u/d quarks.
In the NNPDF2.3qed, the initial neutron's quark/gluon PDFs are related to the proton's with the isospin symmetry, while the violation only occurs at a higher scale as a result of the QED evolution [14].However, the ISV can emerge at the starting scale µ 0 as well.In MRST2004qed, the initial ISV is modeled with the leading-order splitting of the valence quarks since current/constituent masses [13].This idea is inherited in the CT14qed PDFs, while a universal low cutoff scale Q cut = 71 MeV is adopted [8].In this work, we follow the MMHT2015qed/MSHT20qed treatment [22,23] to parameterize the initial ISV as where q (QED) V,p ∝ P (QED) qq ⊗ q V,p .Note that the quark charge e q is implicitly contained in . The ϵ can be fixed through the momentum sum rule of proton and neutron PDFs, as In Fig. 11, we compare the ratios of neutron's valence quarks to their proton isospin partners with the scale from µ = 1.3 GeV up to 1 TeV, which quantifies the size of isospin symmetry violation.The initial ISV effect is about 0.3% (0.7%) at small (large) x, which is less than the one found in the MMHT2015qed PDFs [22], as the elastic component is not included in Eq. (2.38) here.The ISV effect increases with energy scale, which reaches up 1% (6%) at small (large) x when µ = 1 TeV, as a result of the QED DGLAP evolution, consistent with the CT14qed result [8].In such a way, the ISV corrected neutron's valence quarks are obtained through Eq. (2.37), while the sea quarks are still assumed to obey the isospin symmetric assumption, i.e., Eq. (2.19).
The DGLAP evolution.As explicitly demonstrated in Ref. [21], the QED impact on the fitted χ 2 and PDFs in the global analysis is minimal in the CTEQ-TEA framework.Based on this observation, the photon PDF can be obtained with a few iterations, which is enough to give a consistent co-evolved PDF set.
We sketch our procedure as a flowchart shown in Fig. 12.As discussed in the previous subsection, we start with the CT18 proton PDFs and obtain the neutron PDFs with the isospin symmetry approximation.Then, the neutron's photon content at the starting scale xγ(x, µ 2 0 ) can be constructed using the LUXqed formalism.Afterward, the isospin violation induced by the QED effect can be determined correspondingly through Eq. (2.37).The momentum sum rule is enforced by taking out the corresponding photon momentum fraction from the gluon component afterward.Then, the PDFs at a high scale can be self-consistently determined through the DGLAP evolution.With respect to the isospin symmetry approximation (ISA), the DGLAP evolved PDFs can be treated as a small perturbation, which will converge very fast within a few iterations, as depicted in Fig. 12.We take the variation of photon and ISV-corrected valence quark PDFs within the permille level as the criterion to determine the convergence.
In Fig. 13, we present the neutron's inelastic photon, charge-weighted singlet Σ e , as well as g, d, u, Σ PDFs in the CT18qed framework in comparison with the CT18lux ones.Similar to the proton case [21], we see the CT18qed photon PDF gets a few percent enhancement at small x, while agrees quite well in the moderate-x region.The low-x enhancement of CT18qed Proton PDFs Isospin symmetry Eq. (2.19) Neutron PDFs f q,q,g/n (x, µ 2 ) LUXqed photon xγ(x, µ 2 0 ) Initial isospin violation: Momentum sum rule: xγ(x, µ 2 0 ) convergence?no Neutron PDFs f q,q,g,γ/n (x, µ 2 ) yes is induced by the DGLAP evolution, which equivalently integrated out the LO pQCD SFs, larger than the higher-order (NNLO) one employed in the LUXqed formalism, as shown in Fig. 14 left panel.In the extremely large-x region, the CT18lux photon becomes significantly larger than the CT18qed one.It is mainly driven by the larger MS conversion term, induced by a smaller high-Q 2 pQCD SF F 2 than the low-Q 2 no-perturbative one5 , as shown in Fig. 14 right panel.Similar behaviors were found in the proton case already [21].In addition, we also compare the evolved photon with and without isospin symmetry approximation (ISA) at the starting scale.The ISV result includes the enforcement of the momentum sum rule as well.We see that the corresponding difference only shows up at a sub-percent level.At a small x, the ISV photon gets slightly smaller with respect to ISA one, mainly resulting from the momentum sum rule.In comparison, the large-x ISV photon gets a small enhancement induced by the larger charge-weighted singlet.In the right panel of Fig. 13, we present the CT18qed evolved PDFs of g, u, d and singlet Σ, normalized to the CT18lux ones.We see the ISV reduces (enhances) u-and d-quark PDFs at large x, as implied in Fig. 11 already.At small x, the singlet and gluon PDFs get a slight reduction, mainly driven by the enforcement of the momentum sum rule, which turns out to be a small effect, consistent with Eq. (2.34).
Similar to the CT18lux variation, we present various low-Q 2 resources contributing to the neutron's photon PDF uncertainties in the CT18qed framework in Fig. 15.Many features are similar to the CT18lux results in Fig. 9, with some minor differences.For example, the variations induced by Ye et al.'s elastic form factors and R1998 ratios remain at the same level.
The dominant uncertainty comes from the resonance structure functions, by switching the CLAS result to the CB21 fit.In comparison, the relative contribution from the HERMES low-Q 2 continuum region (HERMES and sxQ2) is under control in comparison with the CT18lux case, mainly due to the increase of the perturbative inelastic component through the DGLAP evolution.For this reason, we see the q, g PDF induced uncertainty get enhanced, slightly.In the x → 1 limit, the dominant uncertainty resources are the MHO and TM corrections, propagating from the uncertainties at µ = 1.3 GeV.The HT uncertainty for µ = 100 GeV gets reduced with respect to the CT18lux ones, mainly driven by the small contribution to the initial scale.More specifically, when comparing with the CT18lux photon PDF at µ = 3 GeV as in Fig. 9, we see the CT18qed uncertainties from HERMES, sxQ2, Q2PDF, and MHO get smaller, as a result of its under control at µ 0 = 1.3 GeV.For this reason, we will stay with the µ 0 = 1.3 GeV as our CT18qed default choice for neutron's photon PDF, different from the proton treatment [21].

A comparison with other neutron's photon PDF sets
Similar to the CT18qed proposal [21], we treat these nine low-Q 2 non-perturbative resources as orthogonal eigenvectors for the neutron's Hessian error sets.For Ye et al.'s Elastic form factors, HERMES low-Q 2 continuum SFs, and R1998 ratio, we take the two-directional (±) variations as the corresponding eigenvector sets.For the CB21 resonance SFs, sxQ2 SF variation, Q2PDF matching, MHO, HT, and TM corrections, we symmetrize the difference with respect to the central set to obtain the opposite-directional eigenvector sets.The final combined photon PDF uncertainties can be constructed through where N PDF = 28 refers to the CT18 Hessian eigenvector sets, while n low-Q 2 = 9 corresponds the eight low-Q 2 non-perturbative resources.In Fig. 16, we compare the central values and the error bands for neutron's photon and charge-weighted singlet Σ e = i e 2 i (q 2 i + qi ) in the CT18lux, CT18qed, as well as the MSHT20qed PDFs.In comparison with the CT18lux neutron's Σ e PDF, the CT18qed one gets a 1% reduction around x ∼ 0.6, as a result of the photon's takeaway in the q i → q i γ QED splitting.The relative reduction size is smaller than the proton's case [21] because the dominant valence quark d V,n has a smaller electric charge than the corresponding isospin partners u V,p inside of proton.In comparison with the uncertainty of CT18lux photon at µ = 100 GeV, we see the CT18qed one is slightly smaller, both at small x and large x, due to a small reduction of the non-perturbative impact in the propagation during the DGLAP evolution.For this reason, we advocate the CT18qed as the primary choice, which is used in our study of the phenomenological implications (Sec.3).But overall, we get quite a good agreement for the size of the photon PDF uncertainty between CT18lux and CT18qed.
We also compare our photon and charge-weighted singlet PDFs with the MSHT20qed ones in Fig. 16.We remind that similar to our CT18qed framework, the MSHT20qed initializes the photon PDF with the LUXqed formalism at a low scale, µ 0 = 1 GeV, and evolves the PDF to higher scales [23].Both photon and Σ e PDFs between CT18qed and MSHT20qed agree well within the moderate-x region (10 −3 < x < 10 −1 ).In both small-and large-x regions, the MSHT20qed photon PDF is smaller than the CT18qed ones, similar to the proton case as well [21].MSHT20qed small-x (x ≲ 10 −3 ) photon gets softer, mainly driven by its smaller charge-weighted singlet PDF, as shown in Fig. 16 (right).
In comparison, the MSHT20qed large-x (x ≳ 10 −1 ) photon gets smaller than CT18qed, resulted from various resources.First, MSHT20qed's Σ e at large x is noticeably smaller than the CT18qed one, as shown in Fig. 16 (right).It can be understood in terms of the more flexible parameterization of the MSHT20 fit, which allows a different high-x power for the down quark from the up-quark one, while the CT18 parameterization fixes high-x powers for the up and down quarks to be the same.As a result, the MSHT20qed's neutronto-proton ratio of Σ e approaches to zero while CT18 one does not.Second, as noticed in the CT18qed proton PDF study [21], different treatments of the photon initialization with the LUXqed formalism have been adopted in CT18qed and MSHT20qed.Similar to the MMHT15qed [22], MSHT20qed takes an µ 2 0 rather than µ 2 0 /(1 − z) with Eq. (2.4), while the integration in µ 2 ∈ [µ 2 0 , µ 2 0 /(1 − z)] interval is done with a stationary approximation that which gives a smaller initial photon xγ(x, µ 2 0 ) than the CT18qed one in the x → 1 limit.In Fig. 16, we also compare the relative error bands of photon and charged-weighted singlet PDFs, which are normalized to the corresponding central sets.We see that both PDF sets (CT18 and MSHT20) give a compatible estimation of the photon uncertainty.In the small x region (x < 10 −2 ), the MSHT20qed photon PDF gets a smaller error band, induced by its smaller charged-weighted singlet error band, as shown in Fig. 16 (right).In the region 10 −2 < x < 0.3, the MSHT20qed gives a slightly larger error, induced by its different error estimation [22,23].Following MMHT2015qed [22], MSHT20qed [23] includes variation of elastic form factors, R1998, CLAS/CB resonance SFs, HERMES fits ("Continuum" there), HT ("Renormalon" there) as well as q, g PDFs.Additionally, different from our smooth transition between the HERMES and resonance region, MSHT20qed takes a threshold W 2 cut = 3.5 GeV 2 and varies it within 3 < W 2 cut < 4 GeV 2 to present a separate error set, which gives a larger uncertainty in this region.When x ≳ 0.65, MSHT20qed freezes the photon PDF variation at the starting scale µ 0 , due to the numerical subtlety mentioned above.As a consequence, the photon PDF uncertainty, only driven by the DGLAP evolution, becomes smaller in this region.In comparison, different from the MSHT stationary approximation in Eq. (2.40), we take the complete µ 2 -dependent SFs in the region µ 2 ∈ [µ 2 0 , µ 2 0 /(1 − z)].It allows us to explore the missing higher order (MHO) uncertainty by varying the scale between µ 2 0 and µ 2 0 /(1 − z) in the LUXqed formalism of Eq. (2.4), which gives an additional error set and drives the increasing of CT18qed photon PDF uncertainty when x → 1.
The neutron's photon momentum fractions at a few scales among CT18lux, CT18qed, and MSHT20qed are compared in Table 1.We see the CT18lux and CT18qed results get a quite good agreement, except that the CT18lux result suffers a slightly larger uncertainty at high scale µ = 100 GeV, due to the mentioned TM corrections.The MSHT20 result is smaller, as a result of its smaller elastic photon, as mentioned in Sec.2.1.Meanwhile, we also include the comparison with the first generation of neutron's photon PDF sets, i.e., CT14qed [8], MRST2004qed [13], and NNPDF2.3qed[14], in the Fig. 17.The central PDF of the MRST2004qed set is taken as an average between the two public sets, corresponding to the current-or constitute-quark mass parameterization, respectively, while the error band is their half difference [13].The CT14qed PDF error PDFs are taken as the 0-th and 11-th sets, which correspond to the initial proton's photon momentum fraction ⟨xγ p ⟩(µ 2 0 ) = 0% and 0.11%, respectively, in the 68% confidence level [8].We see that at a low x, all these three first-generation sets give a softer photon, while the MRST2004qed and NNPDF2.3qedgive a larger photon at large x, different from CT14qed.The NNPDF2.3qed gives a significantly larger error band than other PDF sets, due to its different methodology, which fits the available deep inelastic scattering and Drell-Yan data [14].In Fig. 17, we also include the MSHT20qed in this comparison.
Finally, we compare neutron's and proton's photon and charge-weighted singlet PDFs in both CT18qed and MSHT20qed in Fig. 18.We see that the neutron's photon is significantly smaller than the proton's, resulting from its smaller elastic component as well as the suppression of charge-weighted singlet PDFs.In the small-x limit, the neutron's photon approaches the proton's, driven by the unity sea-quark ratio, which is shown explicitly in both the CT18qed and MSHT20qed scenarios.

Phenomenological implications
The neutron's PDFs are essential for many phenomenological studies, especially for the processes involving nucleus scattering.In this section, we take the neutrino-nucleus W -boson production (WBP) and the GeV axion-like particle production in a high-energy muon beam dump experiment to explore the implications of our neutron's photon PDF.

W -boson production in the neutrino-nucleus scattering
The neutrino-nucleus W -boson production (WBP) is an important scattering process for TeV-PeV neutrino detections [57,58].Previous work shows that the cross sections of WBP are up to ≃ 7.5% of the charged current deep-inelastic scattering, which is the dominant process for TeV-PeV neutrino detections by IceCube [59], KM3NeT [60], Baikal-GVD [61], etc.The rapidly increasing data due to more and more running and proposed experiments (e.g., IceCube-Gen2 [62], P-ONE [63], TRIDENT [64], FASERν [65]) necessitate a precision calculation of WBP.Fig. 19 shows the Feynman diagrams of WBP, which is a photon-induced process.Depending on the momentum of the photon, the scattering is in three different kinematic regimes: 1) coherent, in which the photon couples to the nucleus, which remains intact after the scattering, 2) diffractive, in which the photon couples to a nucleon, which remains intact after the scattering, and 3) inelastic, in which the photon couples to a nucleon or quark and the nucleon breaks after the scattering.The inelastic regime has the largest cross section and the largest uncertainty, mainly due to the uncertainties in the photon PDF [66].In Fig. 20, we show the improvement in the precision of the WBP inelastic cross section, thanks to the improvement in the precision of the photon PDF.We follow the computational procedure in Ref. [66].The central factorization and renormalization scales are chosen to be the partonic collision energy √ s νγ , and the corresponding uncertainty is quantified with the 9-point variation by a factor of two.Here we show the WBP cross section of ν e scattered with a neutron target, which is the main focus of this work.The same conclusion can be drawn for proton and for the other lepton flavors.
The upper left panel of Fig. 20 shows the cross sections (divided by neutrino energy to make the plots more visible).The upper right panel shows the cross section ratios with uncertainties and the lower panel shows the relative uncertainties.We see that among the first-generation photon PDF sets, the CT14qed gives the smallest error bands, which was used in Ref. [66].The second-generation photon PDFs (MSHT20qed and CT18qed) give much more precise predictions of WBP inelastic cross section (≲ 10% for E ν ≳ 10 4 GeV as in the right panel).In the E ν ≲ 10 4 GeV regime, the central of MSHT20qed prediction gets smaller than the CT18qed one, due to its smaller large-x photon PDF, as shown in Figs.16-18.With respect to the CT18qed PDF uncertainty, the MSHT20qed one is also smaller and decreases slightly in the small E ν direction, due to its different treatment of photon PDF uncertainty, as discussed in Sec.2.4.
The updated cross sections of WBP using the state-of-the-art CT18qed photon PDF will be uploaded to the GitHub repository under this link .

An Axion-like particle in a muon beam dump experiment
The axion, as a CP-odd scalar particle, was originally postulated by Peccei and Quinn [67] to resolve the strong charge-parity (CP) problem in quantum chromodynamics [68,69].Afterward, many beyond-the-Standard-Model (BSM) theories predict axion or axion-like particles (ALPs) in a wide mass range [70][71][72].See Refs.[2,73] for the latest progress and review for axions and other similar particles.In this subsection, we focus on an ALP (a) with a two-photon interaction with an effective Lagrangian where F µν is the electromagnetic field-strength tensor, with the dual F µν = ϵ µνλρ F λρ /2, and m a is the axion mass.As a consequence, this ALP can be produced through a photon fusion, as shown in Fig. 21, which happens in many high-energy collider environments.For simplicity, we consider a one-side neutron's photon involvement, which can be directly measured at a future high-energy muon beam dump experiment [74] or the muon-ion collider [75].In Fig. 22, we present the ALP production cross sections at a high-energy muon beam dump experiment, with the beam energy E µ = 1.5 TeV, which corresponds to a 3 TeV muon collider [76,77].The corresponding center-of-mass collision energy is where m N is the mass of a nucleon.The calculation is done at the leading order with Mad-Graph [78] interfacing with the UFO [79] model file [80] generated with FeynRules [81].We take f a = 1 TeV as a benchmark and show the ALP production cross section's dependence on the ALP mass from GeV up to m a = 40 GeV, which mainly probes the photon PDF in the large x region.The photon radiated off the neutron target corresponds to the neutron's photon PDF, while the muon's photon is taken from the improved Weizsaecker-Williams approximation [82].The central factorization and renormalization scales are taken as the axion mass, µ F = µ R = m a , while the scale uncertainty is quantified with the 9-point variation.The axion decay width is self-consistently determined solely through the a → γγ channel as We first compare the calculation with only including s-channel diagram and the one with the complete s, t, u ones in the upper left panel of Fig. 22.We see that in the whole m a region, the result with s, t, u channels coincides with the one with the s-channel only very well.This can be understood from the resonance effect, with the negligible contribution from the t, u-channel diagrams.Similar behavior has been observed in the electron-ion collider (EIC) scenario [83].However, some potential additional decay channels, e.g., into invisible dark matter, a → χ χ, can enlarge the axion decay width Γ a , which will signify the destructive interference between the s and t, u channels in Fig. 21 and reduce the production cross section as a consequence.
The PDF and scale uncertainties are presented in Fig. 22, which shows that the scale uncertainty dominates at low m a region, while PDF uncertainty takes over at large m a .The comparison among various neutron's photon PDF results is presented in the rest plots of Fig. 22.We see that with an increase of ALP mass, the ALP production cross section drops very quickly, from about 1 pb down to 10 −4 pb.The scale uncertainty decreases from 60% down to 2%, mainly driven by the decrease of the PDF variation as a result of a smaller high-scale strong coupling.
On the lower left panel, we normalize the cross sections to the CT18qed prediction and compare the corresponding ratios as well as the uncertainties.We see that the MSHT20qed result is in good agreement with the CT18qed result at a low m a when m a ≲ 10 GeV, while gradually gets a significantly smaller when m a ≳ 20 GeV, even up to 60% when m a = 40 GeV, mainly driven by its much smaller photon PDF at large x, as shown in Figs.16-18.As a reference, the CT18qed PDF uncertainty stabilizes at a few percent when m a < 30 GeV, which increases quickly up to 50% when m a = 40 GeV.In comparison, the MSHT20qed uncertainty stays at a 10% level and even decreases a little when ALP mass approaches m a = 40 GeV, which indicates its different large-x photon PDF variation, as discussed in Sec.2.4.
We also include the first generation of photon PDFs in this comparison.We see that the MRST2004qed prediction gets 50% larger than the CT18 benchmark, with a much larger PDF uncertainty (∼ 40%), mainly propagated from the large x PDFs, as shown in Fig. 16.The CT14qed prediction is about 60% ∼ 80% of the CT18qed one, with about 40% variation with respect to its central one.In comparison, both the central prediction and uncertainty band of the NNPDF2.3qedPDF can go much larger than the other PDF sets.

Conclusions
In this work, as a follow-up to our CT18qed photon PDF study [21], we complete it by calculating the neutron's photon content with the latest LUXqed formalism [19,20], which determines the photon PDF with the precisely-measured structure functions.Similar to the proton case, the neutron's PDF is determined in two methodologies, i.e., CT18lux and CT18qed.In the CT18lux approach, the photon PDF is directly calculated at any scale by applying the LUXqed master formula, Eq. (2.4).In comparison, the CT18qed photon PDF is initialized at the starting scale, µ 0 , with the LUXqed, and evolves to a higher scale with the DGLAP evolution, which has reached the NNLO QCD and NLO QED accuracy.
To determine our CT18 neutron's elastic photon PDF, we adopt the neutron's electromagnetic form factors extracted from the global neutron data [3] instead of from the simple Galster parameterization [36,38] used in MMHT2015qed [22] and MHST20qed [23] PDF sets.Moreover, instead of integrating the elastic form factor up to the initial scale µ 0 and evolve to high scale as MMHT2015qed [22] and MHST20qed [23], we take the full LUXqed formalism to determine the elastic component at all scale, which ends up a significant larger elastic photon PDF as shown in Fig. 3.At different scales, the CT18 elastic photon more or less remains the same, with a minor difference coming from the scale evolution of the QED coupling constant.
The inelastic photon PDF, with the CT18lux methodology, is determined using the LUXqed formalism, which integrates the structure functions (SFs) over all the (x, Q 2 ) kinematic regions, as shown in Fig. 1.In the perturbative QCD region with Q 2 > Q 2 PDF , we relate the neutron's quark-gluon PDFs to the proton ones with isospin symmetry approximation.The low-Q 2 structure functions are taken from the experimental measurements and the corresponding fits.In the low-Q 2 continuum region, the HERMES GD-11P and GD-11D fits [29] as well as the R1998 fit [39] of the longitudinal-to-transverse cross-section ratio R L/T = σ L /σ T are adopted.The resonance structure functions are determined through the CLAS [32] or Christy-Bosted (CB) fit [33][34][35] in the resonance region.We take a smooth matching of the SFs in the transition regions between the resonance and low-Q 2 HERMES continuum, as well as the one between the low-Q 2 HERMES and high-Q 2 pQCD continuum.Various low-Q 2 nonperturbative resources, which induce the variation of the inelastic photon are carefully examined, including the CLAS/CB resonance fits, the HERMES, R1998, pQCD SFs in the small-x and small-Q 2 region (sxQ2), higher twist (HT) and target-mass (TM) corrections, PDF matching (Q2PDF), as well as the missing higher order (MHO) uncertainty.It turns out that the resonance SFs (CLAS) can introduce the largest variation at a low scale in the moderate-x range, while the HERMES extrapolation (sxQ2) takes over in the low-x region, both of which dies out at a high scale, while the TM corrections induce a large uncertainty at a large scale and large x.
With the CT18qed methodology, the neutron's inelastic photon is initialized at a low scale µ 0 = 1.3 GeV and evolves to high sales with the NNLO QCD and NLO QED DGLAP equation.With an iterative approach (Fig. 12), the momentum sum rule and the isospin symmetry violation (ISV) are self-consistently determined, which turn out to have a negligible effect on the final photon PDF, as shown in Fig. 13.Similar to the CT18lux, the variation of low-Q 2 nonperturbative resources is examined, which turns out that the CT18qed photon accumulates advantages of the large-scale suppression of the resonance SF effect as well as the small TM corrections to the low-scale PDFs.
The comparison of neutron's photon in the second-generation QED PDFs, i.e., CT18lux, CT18qed, MSHT20qed [23] is performed in Fig. 16.In comparison with the CT18lux, the CT18qed photon PDF gets enhancement at small-x, as a result of the DGLAP evolution, which equivalently integrates out larger LO SFs, with respect to the high-order one in the LUXqed formalism.The larger high-x CT18lux photon is induced by the larger MS matching term induced by a smaller high-Q 2 perturbative F 2 than the low-Q 2 non-perturbative one at the initialization scale in CT18qed.In comparison with the error band, the CT8lux one gets slightly larger error bands, due to a worse control of the TM corrections.For this reason, we advocate the CT18qed as a primary use for phenomenological calculations in related studies.
The MSHT20qed gets a good agreement with the CT18qed in the moderate-x region (10 −3 ≲ x ≲ 10 −1 ), while smaller in both small-and large-x regions.The smaller MSHT20qed small-x photon is due to its smaller charge-weighted singlet, while the smaller large-x photon is induced by both the noticeably smaller charge-weighted singlet and its stationary approximation of Eq. (2.40).The photon PDF uncertainty of MSHT20qed is comparable to CT18qed one, with some moderate difference in the small-x region driven by the charged-weighted singlet while the large-x region driven by the different low-Q 2 non-perturbative SF variations.We also compare our PDFs with the first generation QED PDFs, i.e., MRST2004qed [13], CT14qed [8], and NNPDF2.3qed[14] in Fig. 17, which shows a significant improvement in precision, thanks to the LUXqed formalism.
We also investigate the phenomenological implications of our neutron's photon PDFs with two photon-initiated scattering processes: the neutrino-nucleus W -boson production and the axion-like particle production in a high-energy muon beam-dump experiment.In both scenarios, the cross sections are very sensitive to the neutron's photon PDF at large x.With respect to the first generation of QED PDFs, the second generation ones, both CT18qed and MHST20qed, get a significant theoretical improvement in the PDF uncertainty.In comparison with the CT18qed predictions, the MSHT20qed gets smaller cross sections for the WBP near the energy threshold and for large-mass ALP production at a muon beam dump experiment, due to the corresponding large-x PDF behaviors.
As a companion of this article, the LHAPDF6 [84] grids of the CT18qed PDFs, including the elastic, inelastic components as well as the sum as a total, will be provided through the CTEQ-TEA GitLab repository, https://cteq-tea.gitlab.io/project/00pdfs/.The construction of each eigenvector set is discussed in detail in Sec.2.4.The final PDF uncertainty should be combined in terms of Eq. (2.39).We have re-scaled the n low-Q 2 low-Q 2 nonperturbative error sets with a factor of 1.645 to obtain the 90% CL, respecting the CT Hessian error criterion [56,85,86], for the convenience of combination.Different from the default CT18qed proton QED PDF set which takes the µ 0 = 3 GeV to minimize the nonperturbative PDF uncertainty [21], we take the µ 0 = 1.3 GeV for the neutron one as default, which get better control of the low-Q 2 resource, largely by suppressing the mismatch effect between the HERMES and pQCD continuum regions.
Finally, this and previous work [21] focus on the photon content of free nucleons.The nuclear corrections play a role in various places, such as the neutron's electromagnetic form factors [3] and the resonance structure functions [32][33][34][35].In our treatment, the nuclear uncertainty on the neutron's low-Q 2 SFs in the HERMES continuum region [29] is folded in the HERMES and R1998 error sets.Moreover, the nuclear effects (including shadowing [87,88] or EMC [89]), as well as the electric charge, are expected to have a non-negligible impact on the neutron's inelastic and elastic photon inside a heavy nucleus, such as the lead 208  82 Pb, relevant to the LHC and EIC heavy-ion collisions.This leaves room for future improvement, which can be done with a joint effort with the nuclear PDF fitting group, such as the nCTEQ collaboration.

Figure 1 .
Figure 1.The breakup plane to determine the neutron structure functions, F n 2,L (x, Q 2 ).

Figure 5 .
Figure 5. Upper: the HERMES longitudinal-transverse cross-section sum σ T +L [barn] for the proton (left) and deuteron (right) in the low-Q 2 continuum region.Lower: the derived neutron cross-section sum (right) based on Eq. (2.15), with nuclear corrections folded into the HERMES and R n L/T ratio uncertainties.The shaded area denotes two smooth transition regions.

1 Figure 6 .
Figure 6.The neutron's structure functions (F 2 , F L ) in the smooth transition regionQ 2 lo < Q 2 < Q 2 PDF from the low-Q 2 HERMES to high-Q 2 pQCD continuum.The error bands opening to the low-Q 2 direction are taken from the HERMES data with Eq. (2.16), while the high-Q 2 bands indicate the pQCD PDF uncertainty.

Figure 8 .
Figure 8.The higher-twist (left) and target-mass (right) corrections to the structure functions F 2,L (x, Q 2 ) for a few representative Q 2 in the high-Q 2 continuum region.

Figure 12 .
Figure 12.The flowchart to determine the neutron's inelastic photon with the CT18qed (DGLAP evolution) approach.

Figure 13 .Figure 14 .
Figure13.The PDF ratios of neutron's γ inel , Σ e (left) and g, u, d, Σ (right) at µ = 100 GeV in the CT18qed framework with respect to the CT18lux ones, with isospin symmetry approximation (ISA) or violation (ISV) at starting scale.The "ISV" includes the enforcement of the momentum sum rule as well.

Figure 16 .
Figure 16.The PDF comparison of the central values (upper) and the error bands (lower) for the neutron's photon (left) and charge-weighted singlet (right) at µ = 100 GeV among CT18lux, CT18qed, and MSHT20qed.

Figure 19 .
Figure 19.The t, u-channel diagrams for the W -boson production in the neutrino-nucleus scattering.

Figure 20 .
Figure 20.Cross sections, σ, of ν e -induced WBP in the inelastic regime on a neutron, (−) ν e n → e ∓ W ± .Upper left: σ/E ν from different inelastic neutron photon PDFs.The two curves in the same color represent the uncertainty at a 68% confidence level.Upper right: cross section ratios with uncertainties between CT18qed and others.Lower: Relative uncertainties, δσ/σ.

Figure 21 .
Figure21.The s, t, u-channel diagrams for the axion-like particle production at a muon beam dump experiment.

Figure 22 .
Figure 22.Axion-like particle production cross sections in the muon-neutron scattering at a E µ = 1.5 TeV muon beam dump experiment.The thin and dashed bands indicate the corresponding PDF CL) and scale uncertainties.The upper left panel compares the single s-channel with complete s, t, u-channel calculations, while the rest include all the s, t, u channels.
[56]re 11.Left: ratios of neutron valence quarks d V,n , u V,n to the proton isospin partners at µ = 1.3, 100, 1000 GeV.Right: the same plot with a linear scale to emphasize the large x behavior.The error bands indicate the envelope of CT18 eigenvector sets[56].