Pion Distribution Functions and Kaon Distribution Amplitude and Functions as a Bound System in 1 + 1 Dimensional QCD
Teruo Kuraiorcid
Tokyo, Japan.
DOI: 10.4236/jmp.2026.171004   PDF    HTML   XML   76 Downloads   347 Views  

Abstract

We obtain pion distribution functions and kaon distribution amplitude and functions based on a bound system in 1 + 1 dimensional QCD. Our pion valence quark distribution functions are comparable to xFitter and JAM, so that behavior of our u-quark distribution functions in the range of medium x to close to x = 1 is close to original E615 results. As an asymptotic limit, an exponent of ( 1x ) of our pion distribution functions is 1. In the case of kaon, the ratio of s-quark distribution functions to u-quark distribution functions behaves within 1σ uncertainty region of the latest JAM results. At the same time, the ratio of kaon u-quark distribution functions to that of pion is comparable to N3 results except the range of x0.9 .

Share and Cite:

Kurai, T. (2026) Pion Distribution Functions and Kaon Distribution Amplitude and Functions as a Bound System in 1 + 1 Dimensional QCD. Journal of Modern Physics, 17, 49-81. doi: 10.4236/jmp.2026.171004.

1. Introduction

Recently, pion valence quark and u-quark distribution functions (π v(u)PDF) have become a hot topic in Hadron physics. π v(u)PDF has long history but we can say that people’s recent concern for π v(u)PDF follows two points. One point is asymptotic behavior and the other one is when this asymptotic behavior starts. Concern to asymptotic behavior maintains fairly long after showing the result that asymptotic behavior of E615 results [1] does not include NLO (next leading order) calculation and its reanalysis results [2] include NLO calculation, which are different. The former asymptotic behavior is exponent of ( 1x ) about 1 (contradict to QCD (quantum chromodynamics) prediction), but the latter one is about 2 (or more) (follow to QCD prediction). For theoretical side, the former asymptotic behavior is supported by constituent quark model [3], Nambu-Jona-Lasino model [4], duality argument [5] and recently, Pasquini et al. [6] and Xie et al. [7], while the latter behavior is supported by Dyson Schwinger equation [8]-[12] and recent Lattice QCD calculation [13] [14]. For the second concern, it starts after phenomenological analysis such as JAM [15] and xFitter [16] (JAM and xFitter are working team’s name), which include NLO consideration, show that, in the range of medium x to close to 1, one half of their π vPDF, i.e., π uPDF, behave as close to that of E615 original analysis data, which means close to linear of ( 1x ) when x approaches 1. Note that both xFitter and JAM analysis are based on what form of π vPDF can generate cross-sections that agree with those of considered experiments. However, phenomenological analysis is not sued to argue asymptotic behavior because data set becomes unreliable when x approaches 1 (cross-section itself becomes obscure and number of data is small). In addition, even Dyson Schwinger model [11], which employs inhomogeneous Bethe-Salpeter equation, shows that the actual asymptotic behavior appears at really large Q, that is, x is far closer to 1 than commonly considered.

For kaon distribution functions case, in the absence of empirical information, to separate the quark flavor PDFs in the kaon, first-principle lattice QCD simulations can be used to complement the experimental measurement. For this approach (Lattice QCD), ETM collaboration [17] shows momentum fraction of valence quarks and that of gluons in both pion and kaon. There are only two available experiments for extracting kaon PDFs, i.e., NA3 data [18] and J/ψ production data [19] (these experiments did not show s-quark valence structure). To extract and compare to kaon PDFs from these experiments, we refer, for example, the article of Xu et al. [20] and statistical model by Bourrely et al. [21]. Recently, JAM collaboration [22] shows a first combined QCD analysis of the experimental cross-sections and lattice moments to extract the valence quark PDFs in the pion and kaon and their impact on the gluon momentum fractions. In this paper, in Section 2, we first show π v(u)PDF based on our approximate pion distribution amplitude derived in previous paper [23]. In Section 3, we derive an approximate kaon distribution amplitude and functions which can be compared to that of pion by using nonzero mass solutions for our ‘tHooft model (bound system in 1 + 1 dimensional QCD), same as pion case (zero mass solutions) as shown in [23]. Because s-quark is much heavier than u-quark, we employ the argument of massive quark case as shown in ref. [24], so that we cannot obtain exact solutions in coordinate space. Instead, we set up inhomogeneous integro-differential equations.

2. Pion Distribution Function

We showed a pion distribution amplitude in previous paper [23]. Here, we describe pion distribution functions based on our approximate pion distribution amplitude. Our approximate pion distribution amplitude is described as

| F 3 ( x ) |= 1x x ( 1exp( x 1x ρ )( sin( x 1x ρ )+cos( x 1x ρ ) ) ) (1)

where ρ= u 2| α | , | α |= g 2 16 and x is fraction defined in the region as x[ 0,1 ] , and g 2 is a coupling constant.

The valence quark distribution functions f v ( x ) is defined as

f v ( x ) = normalized distribution amplitude multiplied by x (2)

Thus, in our case, our distribution functions is obtained as

f 3 v ( x )=x | F 3 ( x ) | 0 1 d x | F 3 ( x ) | (3)

Figure 1 shows our pion valence quark distribution functions obtained by Equation (3). Pion valence quark distribution functions with ρ=0.5 and ρ=3 correspond to twice larger than those shown in Figure 2(A) and Figure 2(B) in Raya et al.’s paper [25]. According to Raya et al., a distribution function (π uPDFs) shown in Figure 2(A) represents rest frame case and that shown in Figure 2(B) represents deep inelastic scattering case that is able to compare to E615 data. The reason why our peak values are twice larger than that of common pion u quark distribution functions is that our valence quark distribution functions represent q q ¯ system because they are derived from a charged pion wave function, which represents u d ¯ (or u ¯ d ). Because pion is the Nambu-Goldstone boson mode of QCD, its mass becomes zero at chiral limit (massless quark case). Our wave function is derived in this situation as shown in [23]. Twice larger value is shown both in xFittar [16] and JAM [15] analysis. Important point is that both xFittar and JAM derive a valence quark distribution function by fitting to cross-section data (xFittar used E615, NA10 (286 (GeV) and 194 (GeV)) and WA70 data and JAM used same data as xFittar except Hera data instead of WA70). Thus, in order to compare to common E615 data, we have to multiply a factor 1 2 . Pasquini et al. also show this procedure in ref. [6].

Figure 1. Pion valence quark distribution functions (blue curve ρ=0.5 case and red curve ρ=3 case).

Both xFittar and JAM show that peak value is around 0.7 to 0.8 and peak point is about x = 0.4. Our ρ=3 valence quark distribution functions shows similar results. If we can consider that rest frame distribution functions corresponds to that of initial scale case, we can set 0.85 (GeV) for initial scale Q 0 of our ρ=0.5 distribution functions and that value is adopted by Pasquini et al. [6]. Pasquini et al. describe an evolution pion u-quark distribution functions at Q 2 =27 (GeV2), i.e., Q=5.2 (GeV). Note that Pasquini et al. use notation μ 0 , μ 2 instead of Q. Thus, momentum acceleration is about 6 times as much. We use this consideration to obtain ρ=3 valence quark distribution functions. In ‘tHooft model, mass square is propotional to g 2 /π . This means that dimension of g 2 is (GeV)2, so that dimension of u in ρ is same as momenta (GeV) because ρ is dimensionless and we set c==1 (c is velocity of light and is Plank constant divided by 2π).

Thus, momentum acceleration means that u value is 6 times as much. This dimensional consideration can be confirmed as follows. In previous paper [23], u is set after changing variable as x= e i 3π 4 (or e i π 4 ) 1 | α | z (or z ¯ ) (this x is space coordinate), so that u is dimensionless. In this case, | q | u 2| α | is dimensionless. Then, to obtain x 1x (this x is fraction), | q | is described as (Momentum dimensional value) momentum/ Max momentum 1 momentum/ Max momentum = (momentum dimensional value) q ¯ 1 q ¯ . Then, q ¯ is dimensionless and can be denoted as x (fraction). Now, u is multiplied by this (momentum dimensional value). Thus, dimension of u in ρ is momentum (GeV).

Note that common method to obtain an accelerated distribution function which can be compared to experimental data is adopting DGLAP evolution equation [26] from initial scale distribution function. Wu et al. [27] mention that, in the deep inelastic scattering, the hadron can be regarded as moving with an infinite momentum frame (IMF). In IMF, the hadron is moving with infinite four momentum (E, P) in the z direction. Then, they show that rest frame four-momentum ( p 0 ,p ) transfers to four-momentum ( k 0 ,k ) by inverse Lorenz boost and that rest frame p can be expressed by IMF variables ( ζ= k z P , k ) when taking P , where k z is longitudinal momentum of quarks in IMF (momentum consideration). Using this expression of wave function, they obtain distribution amplitude at initial momentum scale. Then, they obtain distribution function by DGLAP evolution equation (essentially coupling constant consideration). This means that they consider two processes that are momentum consideration and coupling constant consideration to obtain distribution functions from rest frame wave function. In our case, rest frame wave function is obtained in 1 + 1 dimensions. As shown in ref. [23], the proper distribution amplitude in rest frame is obtained by Fourier Transform with choosing ρ value that is choosing coupling constant g 2 value as fixing momentum u (coupling constant consideration). Distribution functions of deep inelastic scattering is expressed as in moving frame so that change of ρ value needs change of momentum u because coupling constant g 2 is already chosen (momentum consideration). Thus, we also consider two processes that are momentum consideration and coupling constant consideration as same as Wu et al.

Figure 2 shows our pion u-quark distribution functions compared to E615 original analysis data [1]. Here, we use the following relation equation between valence quark distribution functions ( q q ¯ system) and u-quark distributions (q only).

Figure 2. Pion u quark distribution functions: pion u-quark (blue curve), E615 data (red dots).

f u ( x )= 1 2 f v ( x ) (4)

In our case,

f 3 u ( x )= 1 2 f 3 v ( x ) (5)

3. Kaon Distribution Function

In order to derive kaon distribution functions, we need to obtain kaon distribution amplitude as shown in the case of that of pion [23]. For pion, as mentioned before, we used the fact that pion mass is zero in the chiral limit (massless quarks) because of pion is the Nambu-Goldstone boson mode. Thus, we used a wave function of our ‘tHooft model with zero mass case. For kaon, although kaon is also considered as the Nambu-Goldstone boson mode of QCD [28], we use the same consideration to construct kaon mass spectrum in 3 + 1 dimension massive quark case [24] instead of using chiral limit consideration. Main reason is following. We cannot express exact next mass besides zero-mass that should be kaon mass [29]. Thus, we prefer to representation of distribution amplitude and functions without explicit dependence of β . In fact, we can obtain distribution amplitude and functions without dependence of β as shown in later by using this consideration. In ref. [24], to obtain kaon mass spectrum, we used a wave function of f 0 ( 500 ) for massless quark case (chiral limit case) and applied the first order perturbation by considering quark mass term as perturbative Hamiltonian. Thus, in this time, we first construct an equation of motion with massive quark in 1 + 1 dimension case. Recalling the fact that distribution amplitude corresponds Fourier Transform of the space coordinate wave function [23], we do not have to consider the corresponding eigenvalue (mass) but need only space coordinate wave function. Thus, we can consider the quark mass term as inhomogeneous part of inhomogeneous second order integro-differential equation.

Before constructing an equation of motion with massive quark case in 1 + 1 dimensions, we check our equation motion in 3 + 1 dimensions by comparing to that of Suura. Suura’s definition of gauge invariant operator q ηξ ( 1,2 ) is described as [30]

q ηξ ( 1,2 )= T r c [ Pexp( ig 1 2 d x A a ( x )( λ a 2 ) ) ] q η ( 1 ) q ξ ( 2 ) (6)

where T r c denotes trace of color spin a and P denotes path ordering with straight line and η,ξ denote Dirac indices. Note that actual Suura’s notation of Dirac indices are α and β .

Note that Suura defined Dirac indices in QED case and Trace for color and path ordering were represented in QCD case. Note that Suura’s representation of path ordering is [ ] + and Trace of color is explained by words.

Our case of gauge invariant operator q ηξ ( 1,2 ) is, for example, described as in ref. [31] as

q ηξ ( 1,2 )= T r c q ξ ( 2 )Pexp( ig 1 2 d x A a ( x )( λ a 2 ) ) q η ( 1 ) (7)

Comparing Equation (6) and Equation (7), only difference is position of q ξ ( 2 ) . If we move this quark field of Equation (7) to the same position of Equation (6), minus sign will appear because of anti-commutation of quark and anti-quark fields. This means that only difference between Equation (6) and Equation (7) is sign. Then, first we compare kinetic terms.

Suura’s definition of kinetic terms is

i α L ( 1 )i α R ( 2 )

where α L( R ) means α operating from left (right) on q( 1,2 ) .

This means that kinetic terms become

kinetic terms=i α ( 1 )q( 1,2 )q( 1,2 )i α ( 2 ) (8)

When we consider relative coordinate as r  =  x ( 2 ) x ( 1 ) , Equation (6) becomes

kinetic terms=[ q s ( 1,2 ),i α ( r ) ]=[ i α ( r ), q s ( 1,2 ) ] =[ i α ( r ), q k ( 1,2 ) ] (9)

where q s ( 1,2 ) and q k ( 1,2 ) denote Suura’s and ours of gauge invariant operator, respectively.

Equation (9) shows that our definition of kinetic terms is also derived from Suura’s one.

For massive quark case, Dirac equation becomes as

i q t =i α β ¯ mq (10)

i q t = q i α + q m β ¯ (11)

Note that normal notation of β ¯ is β .

Equation (10) and Equation (11) are obtained from the representation of Dirac equation as

( γ μ μ +m )q=0 (12)

Adopting metric is η 00 =1 , η 11 = η 22 = η 33 =1 .

Adopting γ -matrices are

γ 0 =( i )( 0 σ 0 σ 0 0 ), γ k =( i )( 0 σ k σ k 0 ), α k = γ 0 γ k and β ¯ =i γ 0

This description is following Weinberg’s one [32].

For Suura’s definition of gauge invariant operator case, massive quark terms become from Equation (10) and Equation (11) as

massive quark terms= m 1 β ¯ q s ( 1,2 )+ q s ( 1,2 ) m 2 β ¯ (13)

Thus, for our definition of gauge invariant operator case, massive quark terms become

massive quark terms= m 1 β ¯ q k ( 1,2 ) q k ( 1,2 ) m 2 β ¯ (14)

Decomposition of gauge invariant operator is

q=1 q 0 +( i α r ^ ) q 1 +β q 2 +β( i α r ^ ) q 3 (15)

Then, massive quark terms become

Unit matrix component term: ( m 1 m 2 ) q 2 (16)

( i α r ^ ) component: ( m 1 + m 2 ) q 3 (17)

β component: ( m 1 m 2 ) q 0 (18)

β( i α r ^ ) component: ( m 1 + m 2 ) q 1 (19)

Thus, by using anti-particle argument mentioned in ref. [24], massive quark terms affecting our equation of motion are following:

( i α r ^ ) component: ( m 1 + m 2 ) q 3 (20)

β( i α r ^ ) component: ( m 1 + m 2 ) q 1 (21)

This is exactly same form obtained in ref. [24].

Above results are important to construct an equation of motion with massive quark in 1 + 1 dimensions because we do not have proper Dirac equation in 1 + 1 (or 2) dimensions. We have to employ Dirac equation with massive quark in 3 + 1 dimensions to that in 1 + 1 dimension as an analogous form.

Thus, our equation of motion with massive quarks in 1 + 1 dimensions becomes

i q ˙ ( 1,2 )=iα( 2 )q( 1,2 )q( 1,2 )iα( 1 )+g 1 2 dx q E ( 1,2:x ) + m 1 q( 1,2 )q( 1,2 ) m 2 (22)

where

q E ( 1,2:x )= T r c q ( 2 )U( 2,x ) E a ( x )U( x,1 )q( 1 ) (23)

U( 2,1 )=Pexp( ig 1 2 dx A a ( x ) λ a 2 ) (24)

Note that E a ( x ) and A a ( x ) indicate E 1 a ( x ) and A 1 a ( x ) , respectively and that Equation (22) except mass terms is same as equation of motion adopting in previous paper [23].

We employ the metric system and γ -matrices in 1 + 1 dimensions as follows.

g 00 =1, g 11 =1

γ 0 =( 0 1 1 0 ), γ 1 =( 0 1 1 0 ),α= γ 0 γ 1 =( 1 0 0 1 )= σ 3

Note that sign of metric system is different from that in ref. [23] because of consideration of massive quark terms as mentioned before, and that we employed γ -matrices defined by Casher et al. [33].

The Bethe-Salpeter like amplitude is defined by sandwiching between vacuum state and physical state as

Χ( 1,2 )=0|q( 1,2 )| Phy (25)

Because decomposition of gauge invariant operator in 1 + 1 dimensions is

q=1 q 0 +i σ 3 q 1 + σ 2 q 2 + σ 1 q 3 (26)

After taking center of mass coordinate and relative coordinate as described in ref. [23], the Bethe-Salpeter like amplitude with massive quarks becomes

P 0 Χ 0 ( r )=i P 1 Χ 1 ( r ) (27)

P 0 Χ 1 ( r )=i P 1 Χ 0 ( r )( m 1 + m 2 ) Χ 3 ( r ) (28)

W 0 Χ 2 ( r )=2 Χ 3 r g 2 8π d x | r x || x | x Χ 3 ( r x ) (29)

W 0 Χ 3 ( r )=2 Χ 2 r + g 2 8π d x | r x || x | x Χ 2 ( r x )( m 1 + m 2 ) Χ 1 ( r ) (30)

Note that we omit O( 1 N ) term because we are considering large N limit case and that when P 1 =0 (rest frame case), P 0 should be W 0 .

From Equation (27) and Equation (28), Χ 1 ( r ) becomes

Χ 1 ( r )= P 0 ( m 1 + m 2 ) P 0 2 P 1 2 Χ 3 ( r ) (31)

Substituting Equation (31) into Equation (30), Equation (30) becomes

W 0 Χ 3 ( r )=2 Χ 2 r + g 2 8π d x | r x || x | x Χ 2 ( r x )+ P 0 ( m 1 + m 2 ) 2 P 0 2 P 1 2 Χ 3 ( r ) (32)

Using Equation (29) and Equation (31) and taking Χ ± ( r )= Χ 3 ( r )±i Χ 2 ( r ) , we obtain following equations:

W 0 Χ ( r )=2i Χ r +i g 2 8π d x | r x || x | x Χ ( r x ) + P 0 ( m 1 + m 2 ) 2 P 0 2 P 1 2 Χ 3 ( r ) (33)

W 0 Χ + ( r )=2i Χ + r i g 2 8π d x | r x || x | x Χ + ( r x ) + P 0 ( m 1 + m 2 ) 2 P 0 2 P 1 2 Χ 3 ( r ) (34)

As described in ref. [23], we solve Equation (33) first. Because it is difficult to solve Equation (33) exactly, we use the following approximation. For the last term, Χ 3 ( r ) , (depending quark mass term), we take mass zero ( W 0 =0 ( | β |=0 )) solution derived in ref. [23] which is related to our pion wave function in 1 + 1 dimensions. Then, we consider the approximated last term of Equation (33) as an inhomogeneous term and after taking derivative with respect to r and some manipulation including Sokhotsky formula to singular integral equation [34] described in ref. [23], we can construct the inhomogeneous second order differential equation as

2 ϕ ( + ) x 2 +[ ( W 0 2i ) g 2 8 1 2i x ] ϕ ( + ) x +[ ( g 2 8π ) g 2 8 1 2i ] ϕ ( + ) =λ c 1 ϕ ( 0 ) ( + ) x (35)

2 ϕ ( ) x 2 +[ ( W 0 2i )+ g 2 8 1 2i x ] ϕ ( ) x +[ ( g 2 8π )+ g 2 8 1 2i ] ϕ ( ) =λ c 2 ϕ ( 0 ) ( ) x (36)

where λ= 1 2i P 0 ( m 1 + m 2 ) 2 P 0 2 P 1 2 , c 1 and c 2 are constants determined later.

Note that we change the notation of space coordinate r to x in Equation (35) and Equation (36) because we are working in the framework of 1 + 1 dimensions. From now on, we use this notation for space coordinate.

Note that Equation (35) and Equation (36) are composed by the quantities when x asymptotically approaches real axis in the upper-half hemisphere and in the lower-half hemisphere, respectively.

Homogeneous parts of equations for Equation (35) and Equation (36) are obtained by setting left-hand side of equations be zero for both cases as

2 ϕ ( + )( 0 ) x 2 +[ ( W 0 2i ) g 2 8 1 2i x ] ϕ ( + )( 0 ) x +[ ( g 2 8π ) g 2 8 1 2i ] ϕ ( + )( 0 ) =0 (37)

2 ϕ ( )( 0 ) x 2 +[ ( W 0 2i )+ g 2 8 1 2i x ] ϕ ( )( 0 ) x +[ ( g 2 8π )+ g 2 8 1 2i ] ϕ ( )( 0 ) =0 (38)

One of the solutions is given in ref. [23] for both equations. These solutions are as follows:

ϕ ( 1 ) ( + )( 0 ) ( x )= e 1 4 α x 2 e 1 2 βx 2 1 4 i π ( α x β α ) 1 2 W 1 4 i π , 1 4 ( 1 2 ( α x β α ) 2 ) (39)

ϕ ( 1 ) ( )( 0 ) ( x )= e 1 4 α x 2 e 1 2 βx 2 1 4 + i π ( i( α x+ β α ) ) 1 2 W 1 4 i π , 1 4 ( 1 2 ( α x+ β α ) 2 ) (40)

where α= g 2 8 1 2i , β= W 0 2i .

As shown in Appendix C of ref. [23], the characteristic part of solution described in Equation (39) is Weber function D λ ¯ ( t ) type and that in Equation (40) is Weber function D λ ¯ 1 ( it ) type.

Thus, the other solution for Equation (35) must be Weber function D λ ¯ 1 ( it ) type and that for Equation (36) must be Weber function D λ ¯ ( t ) type. Because of λ ¯ =1 2i π for Equation (35) and λ ¯ = 2i π for Equation (36), the other solutions for Equation (35) and Equation (36) become

ϕ ( 2 ) ( + )( 0 ) ( x )= e 1 4 α x 2 e 1 2 βx 2 1 4 i π ( i( α x β α ) ) 1 2 W 1 4 + i π , 1 4 ( 1 2 ( α x β α ) 2 ) (41)

ϕ ( 2 ) ( )( 0 ) ( x )= e 1 4 α x 2 e 1 2 βx 2 1 4 + i π ( α x+ β α ) 1 2 W 1 4 i π , 1 4 ( 1 2 ( α x+ β α ) 2 ) (42)

Then, Equation (39) and Equation (41) are considered as basic solutions for homogeneous second order differential equation described in Equation (35) and Equation (40) and Equation (42) are considered as basic solutions for homogeneous second order differential equation described in Equation (36).

In order to construct particular solutions, we follow the argument of Ince [35]. Our case’s Wronskian are defined as

Wron s ( 1 ) =| ϕ ( 1 ) ( + )( 0 ) ( x ) ϕ ( 1 ) ( + )( 0 ) x ϕ ( 2 ) ( + )( 0 ) ( x ) ϕ ( 2 ) ( + )( 0 ) x | (43)

Wro s ( 2 ) =| ϕ ( 1 ) ( )( 0 ) ( x ) ϕ ( 1 ) ( )( 0 ) x ϕ ( 2 ) ( )( 0 ) ( x ) ϕ ( 2 ) ( )( 0 ) x | (44)

Equation (43) is Wronskian of Equation (35) and Equation (44) is Wronskian of Equation (36).

Using the following relation equation for derivative of Whittaker function W κ,μ ( z ) [36]

z W κ,μ ( z )=( z 2 κ ) W κ,μ ( z ) W κ+1,μ ( z ) (45)

the form of Wron s ( 1 ) and Wron s ( 2 ) are described as

Wron s ( 1 ) = e i| β |x [ W κ 1 ( 1 ) , 1 4 ( i 1 2 ( | α | x | β | | α | ) 2 ) W κ 2 ( 1 ) , 1 4 ( i 1 2 ( | α | x | β | | α | ) 2 ) ×( 2 | α | ( κ 1 ( 1 ) κ 2 ( 1 ) ) ( | α | x | β | | α | ) 2 +i | α | ) +2 | α | ( | α | x | β | | α | ) 2 ( W κ 2 ( 1 ) , 1 4 ( i 1 2 ( | α | x | β | | α | ) 2 ) × W κ 1 ( 1 ) +1, 1 4 ( i 1 2 ( | α | x | β | | α | ) 2 ) W κ 1 ( 1 ) , 1 4 ( i 1 2 ( | α | x | β | | α | ) 2 ) W κ 2 ( 1 ) +1, 1 4 ( i 1 2 ( | α | x | β | | α | ) 2 ) ) ] (46)

where

α=i| α |=i g 2 16 ,β=i| β |=i W 0 2

κ 1 ( 1 ) = 1 4 i π , κ 2 ( 1 ) = 1 4 + i π

Wron s ( 2 ) = e i| β |x [ W κ 1 ( 2 ) , 1 4 ( i 1 2 ( | α | x+ | β | | α | ) 2 ) W κ 2 ( 2 ) , 1 4 ( i 1 2 ( | α | x+ | β | | α | ) 2 ) ×( 2 | α | ( κ 1 ( 2 ) κ 2 ( 2 ) ) ( | α | x+ | β | | α | ) 2 +i | α | ) +2 | α | ( | α | x+ | β | | α | ) 2 ( W κ 2 ( 2 ) , 1 4 ( i 1 2 ( | α | x+ | β | | α | ) 2 ) × W κ 1 ( 2 ) +1, 1 4 ( i 1 2 ( | α | x+ | β | | α | ) 2 ) W κ 1 ( 2 ) , 1 4 ( i 1 2 ( | α | x+ | β | | α | ) 2 ) W κ 2 ( 2 ) +1, 1 4 ( i 1 2 ( | α | x+ | β | | α | ) 2 ) ) ] (47)

where

κ 1 ( 2 ) = 1 4 i π , κ 2 ( 2 ) = 1 4 + i π

According to Ince [35], when the basic solutions of homogeneous part are given as M 1 ( x ) and M 2 ( x ) and inhomogeneous part is given as f( x ) , the particular solution is defined as

Particular solution= x dt f( t ) M 1 ( t ) Wroskian( t ) M 2 ( x ) x dt f( t ) M 2 ( t ) Wronskian( t ) M 1 ( x )

Lower limit of integral is chosen by boundary condition consideration.

Thus, our case of particular solutions becomes

part .sol ( 1 ) = 0 x dt λ c 1 e i 1 2 α t 2 e i 1 2 β t ( | α | t ) 1 2 W 1 4 i π , 1 4 ( i 1 2 | α | t 2 ) Wron s 1 ( 1 ) ( t ) ϕ ( 2 ) ( + )( 0 ) ( x ) 0 x dt λ c 1 e i 1 2 β t ( | α | t ) 1 2 W 1 4 i π , 1 4 ( i 1 2 | α | t 2 ) Wron s 2 ( 1 ) ( t ) ϕ ( 1 ) ( + )( 0 ) ( x ) (48)

part .sol ( 2 ) = 0 x dt λ c 2 e i 1 2 α t 2 e i 1 2 β t ( | α | t ) 1 2 W 1 4 + i π , 1 4 ( i 1 2 | α | t 2 ) Wron s 1 ( 2 ) ( t ) ϕ ( 2 ) ( )( 0 ) ( x ) 0 x dt λ c 2 e i 1 2 β t ( | α | t ) 1 2 W 1 4 + i π , 1 4 ( i 1 2 | α | t 2 ) Wron s 2 ( 2 ) ( t ) ϕ ( 1 ) ( )( 0 ) ( x ) (49)

Note that we omit coefficients of basic solutions in both Equation (48) and Equation (49) because these are inefficient for particular solutions.

Each Wronskian is described as

Wron s 1 ( 1 ) = W 1 4 + 1 π , 1 4 ( i 1 2 ( | α | t | β | | α | ) 2 )( 2 | α | ( 1 2 + 2i π ) ( | α | t | β | | α | ) 3 2 +i | α | ( | α | t | β | | α | ) 1 2 ) +2 | α | ( | α | t | β | | α | ) 3 2 ( W 1 4 + 1 π , 1 4 ( i 1 2 ( | α | t | β | | α | ) 2 ) W 3 4 i π , 1 4 ( i 1 2 ( | α | t | β | | α | ) 2 ) W 1 4 i π , 1 4 ( i 1 2 ( | α | t | β | | α | ) 2 ) W 5 4 + 1 π , 1 4 ( i 1 2 ( | α | t | β | | α | ) 2 ) ) (50)

Wron s 2 ( 1 ) = W 1 4 1 π , 1 4 ( i 1 2 ( | α | t | β | | α | ) 2 )( 2 | α | ( 1 2 + 2i π ) ( | α | t | β | | α | ) 3 2 +i | α | ( | α | t | β | | α | ) 1 2 ) +2 | α | ( | α | t | β | | α | ) 3 2 ( W 1 4 1 π , 1 4 ( i 1 2 ( | α | t | β | | α | ) 2 ) W 5 4 + i π , 1 4 ( i 1 2 ( | α | t | β | | α | ) 2 ) W 1 4 + i π , 1 4 ( i 1 2 ( | α | t | β | | α | ) 2 ) W 3 4 1 π , 1 4 ( i 1 2 ( | α | t | β | | α | ) 2 ) ) (51)

Wron s 1 ( 2 ) = W 1 4 1 π , 1 4 ( i 1 2 ( | α | t+ | β | | α | ) 2 )( 2 | α | ( 1 2 + 2i π ) ( | α | t+ | β | | α | ) 3 2 +i | α | ( | α | t+ | β | | α | ) 1 2 ) +2 | α | ( | α | t | β | | α | ) 3 2 ( W 1 4 1 π , 1 4 ( i 1 2 ( | α | t+ | β | | α | ) 2 ) W 3 4 + i π , 1 4 ( i 1 2 ( | α | t+ | β | | α | ) 2 ) W 1 4 + i π , 1 4 ( i 1 2 ( | α | t+ | β | | α | ) 2 ) W 5 4 1 π , 1 4 ( i 1 2 ( | α | t+ | β | | α | ) 2 ) ) (52)

Wron s 2 ( 2 ) = W 1 4 + 1 π , 1 4 ( i 1 2 ( | α | t+ | β | | α | ) 2 )( 2 | α | ( 1 2 + 2i π ) ( | α | t+ | β | | α | ) 3 2 +i | α | ( | α | t+ | β | | α | ) 1 2 ) +2 | α | ( | α | t+ | β | | α | ) 3 2 ( W 1 4 + 1 π , 1 4 ( i 1 2 ( | α | t+ | β | | α | ) 2 ) W 5 4 i π , 1 4 ( i 1 2 ( | α | t+ | β | | α | ) 2 ) W 1 4 1 π , 1 4 ( i 1 2 ( | α | t+ | β | | α | ) 2 ) W 3 4 + 1 π , 1 4 ( i 1 2 ( | α | t+ | β | | α | ) 2 ) ) (53)

To construct distribution amplitude, that is a function in momentum space, we consider Fourier Transform of particular solution F p ( q ) as shown in ref. [23].

F p ( 1 ) ( q )= dx e iqx 0 x dt c 1 e i 1 2 α t 2 e i 1 2 β t ( | α | t ) 1 2 W 1 4 i π , 1 4 ( i 1 2 | α | t 2 ) Wron s 1 ( 1 ) ( t ) × e i 1 4 | α | x 2 e i 1 2 | β |x ( | α | x | β | | α | ) 1 2 W 1 4 + i π , 1 4 ( i 1 2 ( | α | x | β | | α | ) 2 ) + dx 0 x dt c 1 e i 1 2 β t ( | α | t ) 1 2 W 1 4 i π , 1 4 ( i 1 2 | α | t 2 ) Wron s 2 ( 1 ) ( t ) × e i 1 4 | α | x 2 e i 1 2 | β |x ( | α | x | β | | α | ) 1 2 W 1 4 i π , 1 4 ( i 1 2 ( | α | x | β | | α | ) 2 ) (54)

F p ( 2 ) ( q )= dx e iqx 0 x dt c 2 e i 1 2 β t ( | α | t ) 1 2 W 1 4 + i π , 1 4 ( i 1 2 | α | t 2 ) Wron s 1 ( 2 ) ( t )

× e i 1 4 | α | x 2 e i 1 2 | β |x ( | α | x | β | | α | ) 1 2 W 1 4 + i π , 1 4 ( i 1 2 ( | α | x+ | β | | α | ) 2 ) + dx 0 x dt c 2 e i 1 2 | α | t 2 e i 1 2 β t ( | α | t ) 1 2 W 1 4 + i π , 1 4 ( i 1 2 | α | t 2 ) Wron s 2 ( 2 ) ( t ) × e i 1 4 | α | x 2 e i 1 2 | β |x ( | α | x | β | | α | ) 1 2 W 1 4 i π , 1 4 ( i 1 2 ( | α | x+ | β | | α | ) 2 ) (55)

First, we consider Equation (54).

We work out integral as dx = 0 dx + 0 dx .

For the first integral, changing variables x= x ¯ and t= t ¯ , integral part, I 1 ( 1 ) = 0 dx , becomes

I 1 ( 1 ) = 0 d x ¯ e i| q | x ¯ 0 x ¯ d t ¯ c 1 e i 1 2 α t ¯ 2 e i 1 2 β t ¯ i ( | α | t ¯ ) 1 2 W 1 4 i π , 1 4 ( i 1 2 | α | t ¯ 2 ) Wron s 1 ( 1 ) ( t ¯ ) × e i 1 4 | α | x ¯ 2 e i 1 2 | β | x ¯ i ( | α | x ¯ + | β | | α | ) 1 2 W 1 4 + i π , 1 4 ( i 1 2 ( | α | x+ | β | | α | ) 2 ) = 0 d x ¯ e i| q | x ¯ 0 x ¯ d t ¯ c 1 e i 1 2 α t ¯ 2 e i 1 2 β t ¯ ( | α | t ¯ ) 1 2 W 1 4 i π , 1 4 ( i 1 2 | α | t ¯ 2 ) Wron s 1 ( 1 ) ( t ¯ ) × e i 1 4 | α | x ¯ 2 e i 1 2 | β | x ¯ ( | α | x ¯ + | β | | α | ) 1 2 W 1 4 + i π , 1 4 ( i 1 2 ( | α | x+ | β | | α | ) 2 ) (56)

Note that we take q=| q | in this region by using the argument in ref. [23] and taking 1= e iπ here.

For Wron s 1 ( 1 ) ( t ¯ ) ,

Wron s 1 ( 1 ) ( t ¯ ) = W 1 4 + 1 π , 1 4 ( i 1 2 ( | α | t ¯ + | β | | α | ) 2 )( 2 | α | ( 1 2 + 2i π )( i ) ( | α | t ¯ + | β | | α | ) 3 2 + i | α | ( i ) ( | α | t ¯ + | β | | α | ) 1 2 )+2 | α | ( i ) ( | α | t ¯ + | β | | α | ) 3 2 ×( W 1 4 + 1 π , 1 4 ( i 1 2 ( | α | t ¯ + | β | | α | ) 2 ) W 3 4 i π , 1 4 ( i 1 2 ( | α | t ¯ + | β | | α | ) 2 ) W 1 4 i π , 1 4 ( i 1 2 ( | α | t ¯ + | β | | α | ) 2 ) (57)

W 5 4 + 1 π , 1 4 ( i 1 2 ( | α | t ¯ + | β | | α | ) 2 ) )

Note that ( 1 ) 3 2 = e i 3 2 π =i because we take 1= e iπ in Equation (56).

Equation (57) shows change of variables generates the factor i for Equation (56).

For large | q | case, it is sufficient to consider very small x . Then, integral part of particular solution can be described by series expansion as

I p ( 1 ) = I p ( 1 ) ( x=0 )+ I p ( 1 ) x x+ 1 2 2 I p ( 1 ) x 2 x 2 (58)

Recalling the fact that W 0 is not zero but small ( | β | is small) and that we are considering very small x case, integrand of the first part of particular solutions Integran d 1 ( 1 ) can be expressed as

(59)

Recalling the fact that W κ,μ ( z ) can be expressed by the first term of M κ,μ ( z ) of which exponent of z is dependent of only μ when z is very small ( | α | x+ | β | | α | is very small), we can set W 3 4 i π , 1 4 ( i 1 2 ( | α | x+ | β | | α | ) 2 ) W 1 4 i π , 1 4 ( i 1 2 ( | α | x+ | β | | α | ) 2 ) = Γ( 1+ i π ) Γ( i π ) = i π , W 5 4 + 1 π , 1 4 ( i 1 2 ( | α | x+ | β | | α | ) 2 ) W 1 4 + 1 π , 1 4 ( i 1 2 ( | α | x+ | β | | α | ) 2 ) = 1 1 2 1 π , so that we obtain the second line.

Using the same argument for I 2 ( 1 ) = 0 dx part, we obtain its integrand as

Integran d 2 ( 1 ) ( t ) c 1 e i 1 4 α x 2 ( | α | x ) 1 2 W 1 4 i π , 1 4 ( i 1 2 | α | x 2 ) ( | α | x | β | | α | ) 3 2 W 1 4 + 1 π , 1 4 ( i 1 2 ( | α | x | β | | α | ) 2 )( 2 | α | ( 1 2 i π 1 1 2 i π ) ) (60)

From the integral region, I p ( 1 ) ( 0 )=0 . The particular solution which is composed by multiplying the second term of expansion I p( 2 ) ( 1 ) = I p ( 1 ) x x to basic solution becomes

I p( 2 )( t ) ( 1 ) = c 1 e i 1 4 α x 2 ( | α | x ) 1 2 W 1 4 i π , 1 4 ( i 1 2 | α | x 2 ) ( | α | x+ | β | | α | ) 3 2 W 1 4 + 1 π , 1 4 ( i 1 2 ( | α | x+ | β | | α | ) 2 )( 2 | α | ( i )( 1 2 + 2i π )+i | α | ( i ) ( | α | x+ | β | | α | ) 2 +2 | α | ( i )( i π 1 1 2 i π ) ) x ( | α | x ¯ + | β | | α | ) 1 2 W 1 4 + i π , 1 4 ( i 1 2 ( | α | x+ | β | | α | ) 2 ) =i c 1 e i 1 4 α x 2 ( | α | x ) 1 2 W 1 4 i π , 1 4 ( i 1 2 | α | x 2 )( | α | x+ | β | | α | )x ( 2 | α | ( 1 2 i π 1 1 2 i π ) ) (61)

I p( 2 )( t ) ( 1 ) = c 1 e i 1 4 α x 2 ( | α | x ) 1 2 W 1 4 i π , 1 4 ( i 1 2 | α | x 2 )( | α | x | β | | α | )x ( 2 | α | ( 1 2 i π 1 1 2 i π ) ) (62)

Recalling that we are considering the case that | β | is very small, we obtain the first part of Fourier Transform of the particular solution as

F p( 1 ) ( 1 ) =( 1+i ) 0 dx e i| q |x c 1 e i 1 4 α x 2 ( | α | x ) 1 2 W 1 4 i π , 1 4 ( i 1 2 | α | x 2 ) | α | x 2 ( 2 | α | ( 1 2 i π 1 1 2 i π ) ) (63)

We employ same consideration for the Fourier Transform of the second part of particular solution. Then, F p( 2 ) ( 1 ) becomes as

F p( 2 ) ( 1 ) =( 1+i ) 0 dx e i| q |x c 1 e i 1 4 α x 2 ( | α | x ) 1 2 W 1 4 i π , 1 4 ( i 1 2 | α | x 2 ) | α | x 2 ( 2 | α | ( 1 2 3i π + 1 1 2 i π ) ) (64)

Note that ( 1+i ) term appears in Equation (64) because we take 1= e iπ for Wron s 2 ( 1 ) ( t ¯ ) case as same as that for Wron s 1 ( 1 ) ( t ¯ ) case.

For Equation (55), we employ the same argument used for Equation (54).

Then, after changing variables as x= x ¯ , and t= t ¯ , Wronskian of the first part of the particular solution becomes

Wron s 1 ( 2 ) ( t ¯ ) = W 1 4 1 π , 1 4 ( i 1 2 ( | α | t ¯ | β | | α | ) 2 )( 2 | α | ( 1 2 + 2i π )i ( | α | t ¯ | β | | α | ) 3 2 + i | α | i ( | α | t ¯ | β | | α | ) 1 2 )+2 | α | i ( | α | t ¯ | β | | α | ) 3 2 ×( W 1 4 1 π , 1 4 ( i 1 2 ( | α | t ¯ | β | | α | ) 2 ) W 3 4 + i π , 1 4 ( i 1 2 ( | α | t ¯ | β | | α | ) 2 ) W 1 4 + i π , 1 4 ( i 1 2 ( | α | t ¯ | β | | α | ) 2 ) W 5 4 1 π , 1 4 ( i 1 2 ( | α | t ¯ | β | | α | ) 2 ) ) (65)

Equation (65) shows Wronskian part generates 1 i for integrand for Equation (55) case because we take 1= e iπ instead of e iπ here.

Using the series expansion for integral part as shown in Equation (58), Fourier Transform of the first part of the particular solution of Equation (55) becomes

F p( 1 ) ( 2 ) =( 1i ) 0 dx e i| q |x c 2 e i 1 4 α x 2 ( | α | x ) 1 2 W 1 4 + i π , 1 4 ( i 1 2 | α | x 2 ) | α | x 2 ( 2 | α | ( 1 2 + 3i π 1 1 2 + i π ) ) (66)

F p( 2 ) ( 2 ) =( 1i ) 0 dx e i| q |x c 2 e i 1 4 α x 2 ( | α | x ) 1 2 W 1 4 + i π , 1 4 ( i 1 2 | α | x 2 ) | α | x 2 ( 2 | α | ( 1 2 + i π + 1 1 2 + i π ) ) (67)

To solve Equation (34), we can use same argument for solving Equation (33) as shown before. As shown in ref. [23], basic solutions of Χ + ( x ) can obtain by replacing W 0 to W 0 , that is replacing β to β , in those of Χ ( x ) . That means that we can obtain basic solutions of Χ + ( x ) by replacing β to β in Equation (39) to Equation (42). Because inhomogeneous part is the same as that of Equation (33), all argument for solving Equation (33) can be used by just replacing β to β . Because Fourier Transform of the particular solution for Equation (33) is independent of β shown in Equations (63) - (64) and Equations (66) - (67), Fourier Transform of the particular solution of Equation (34) is composed by the same form as in Equations (63) - (64) and Equations (66) - (67). Recalling the relation Χ 3 = Χ + Χ + 2 , we can construct Fourier Transform of Χ 3 ( x ) as follows.

To invoke the condition that inhomogeneous part is pion wave function, we set c 1 and c 2 as following.

c 1 ( 1 2( 1 2 i π 1 1 2 i π ) + 1 2( 1 2 3i π + 1 1 2 i π ) )= 2 i π Γ( 1+ i π ) e i π 8 (68)

c 2 ( 1 2( 1 2 + 3i π 1 1 2 + i π ) + 1 2( 1 2 + i π + 1 1 2 + i π ) )= 2 i π Γ( 1 i π ) e i π 8 (69)

Note that | α | term is cancelled out.

Then, Fourier Transform of Χ 3 ( x ) can be expressed as

F 3 ( | q | )=( 1+i ) 0 dx e | q |x Γ( 1+ i π ) e i π 8 e i 1 4 α x 2 ( | α | x ) 1 2 W 1 4 i π , 1 4 ( i 1 2 | α | x 2 ) x 2 +( 1i ) 0 dx e | q |x Γ( 1 i π ) e i π 8 e i 1 4 α x 2 ( | α | x ) 1 2 W 1 4 + i π , 1 4 ( i 1 2 | α | x 2 ) x 2 (70)

Note that λ defined in Equations (35) - (36) is omitted because we are interesting in constructing distribution amplitude that is defined by normalized form and that, apart from constant, Equation (70) is exactly same form as Fourier Transform of pion wave function in large | q | case after taking α=| α | e i π 2 in ref. [23] except both integrand are multiplied by x 2 factor. Also note that constant factor e i π 8 and e i π 8 are required to the condition that, for small | q | case, the order of | q | of the first term of F 3 ( | q | ) is linear of | q | (exponent of | q | is 1).

Before using change of variables shown in ref. [23], we employ the representation of W κ,μ ( z ) as following form [36] because we are interested in very large | q | case in which only very small x is sufficient to consider for integral.

W κ,μ ( z )= Γ( 2μ ) Γ( 1 2 μκ ) M κ,μ ( z )+ Γ( 2μ ) Γ( 1 2 +μκ ) M κ,μ ( z ) (71)

where

M κ,μ ( z )= z μ+ 1 2 e z 2 F 1 1 ( μκ+ 1 2 ,2μ+1;z )

F 1 1 is confluent hyper geometric series defined as

F 1 1 ( ζ,ξ;z )=1+ n=1 ζ( ζ+1 )( ζ+n1 ) ξ( ξ+1 )( ξ+n1 ) z n n!

For very small x case, we need only the first term of Equation (71).

Recalling that μ= 1 4 , κ= 1 4 i π for the first integral and κ= 1 4 + i π for the second integral, respectively, Equation (70) becomes

F 3 ( | q | )= 0 dx e i| q |x x 2 F 1 1 ( 1 2 + i π , 1 2 ;i 1 2 | α | x 2 ) 0 dx e i| q |x x 2 F 1 1 ( 1 2 i π , 1 2 ;i 1 2 | α | x 2 ) (72)

Note that ( i 1 2 | α | x 2 ) 1 4 + 1 2 = e i π 8 ( 1 2 | α | x ) 1 2 , ( i 1 2 | α | x 2 ) 1 4 + 1 2 = e i π 8 ( 1 2 | α | x ) 1 2 and those e i π 8 and e i π 8 are cancelled out by our previous setting as mentioned before and that factor ( 1 2 ) 1 2 is omitted because we take normalization later as mentioned before.

For Equation (72), we apply change of integral variables as same as those in ref. [23], which are x= e i 3 4 π z ¯ | α | for the first integral in Equation (72) and

x= e i 1 4 π z ¯ | α | for the second integral in Equation (72), respectively. Using contour integral as shown in ref. [23] and considering the first term of confluent hyper geometric series F 1 1 because we are considering very large | q | case (very small z ¯ case), Equation (72) becomes as

F 3 ( | q | )=( 1+i ) e i 1 4 π | α | 5 2 0 d z ¯ exp( | q | 2| α | z ¯ )exp( i | q | 2| α | z ¯ ) z ¯ 2 ( 1i ) e i 3 4 π | α | 5 2 0 d z ¯ exp( | q | 2| α | z ¯ )exp( i | q | 2| α | z ¯ ) z ¯ 2 (73)

Note that e i 1 4 π in the first line comes from the fact that e i 9 4 π = e i2π e i 1 4 π = e i 1 4 π .

To obtain the final form of approximate form of F 3 ( | q | ) for very large | q | case, we divided integral as

0 dx = 0 u dx + u dx (74)

following the argument in ref. [23]. For very large | q | case, we need only the first term of Equation (74) as shown in ref. [23]. Then, integral part of Equation (73) becomes

I 3 ( 1 ) = 0 u d z ¯ exp( | q | 2| α | z ¯ )exp( i | q | 2| α | z ¯ ) z ¯ 2 (75)

I 3 ( 2 ) = 0 u d z ¯ exp( | q | 2| α | z ¯ )exp( i | q | 2| α | z ¯ ) z ¯ 2 (76)

Denoting μ = | q | 2| α | ( 1i ) and μ + = | q | 2| α | ( 1+i ) , Equation (75) and Equation (76) become

I 3 ( 1 ) = 1 μ e μ u u 2 2 μ 2 e μ u u+ 2 μ 3 ( 1 e μ u ) (77)

I 3 ( 2 ) = 1 μ + e μ + u u 2 2 μ 2 e μ + u u+ 2 μ 3 ( 1 e μ + u ) (78)

Because we are considering in the case of very large | q | , it is sufficient to consider the last terms of both Equation (77) and Equation (78) only. Then, we obtain the final form of approximate form of F 3 ( | q | ) for very large | q | case as

F 3 ( | q | )=( 1+i ) e i 1 4 π | q | 3 1 ( 1i 2 ) 3 1 | α | ( 1exp( | q | 2| α | u )exp( i | q | 2| α | u ) ) ( 1i ) e i 3 4 π | q | 3 1 ( 1+i 2 ) 3 1 | α | ( 1exp( | q | 2| α | u )exp( i | q | 2| α | u ) ) = e i 1 2 π [ 2 e i 1 4 π e iπ 1 | q | 3 | α | ( 1exp( | q | 2| α | u )exp( i | q | 2| α | u ) ) 2 e i 1 4 π e iπ 1 | q | 3 | α | ( 1exp( | q | 2| α | u )exp( i | q | 2| α | u ) ) ] =2 2 1 | q | 3 | α | [ sin( π 4 )exp( | q | 2| α | u )sin( | q | 2| α | + π 4 ) ] =2 1 | q | 3 | α | [ 1exp( | q | 2| α | u )( sin( | q | 2| α | )+cos( | q | 2| α | ) ) ] (79)

For distribution amplitude, we consider absolute value of F 3 ( | q | ) as shown in ref. [23]. Because distribution amplitude is defined by normalized form, the final form of approximate form of | F 3 ( | q | ) | for very large | q | case becomes as

| F 3 ( | q | ) |= 1 | q | 3 ( 1exp( | q | 2| α | u )( sin( | q | 2| α | u )+cos( | q | 2| α | u ) ) ) (80)

For very small | q | case, we have to consider the case that x is sufficiently large so that we first perform integration of integral part of particular solutions. For upper half hemisphere group, integral part of particular solutions, I 1 ( 1 ) and I 2 ( 1 ) , are expressed as

I 1 ( 1 ) = e i | β | 2 x 0 x d t ¯ c 1 e i 1 2 α t ¯ 2 e i 1 2 β t ¯ ( | α | t ¯ ) 1 2 W 1 4 i π , 1 4 ( i 1 2 | α | t ¯ 2 ) Wron s 1 ( 1 ) ( t ¯ ) + e i | β | 2 x 0 x d t ¯ c 1 e i 1 2 α t ¯ 2 e i 1 2 β t ¯ ( | α | t ¯ ) 1 2 W 1 4 i π , 1 4 ( i 1 2 | α | t ¯ 2 ) Wron s 1 ( 1 ) ( t ¯ ) (81)

I 2 ( 1 ) = e i | β | 2 x 0 x d t ¯ c 1 e i 1 2 β t ¯ ( | α | t ¯ ) 1 2 W 1 4 i π , 1 4 ( i 1 2 | α | t ¯ 2 ) Wron s 2 ( 1 ) ( t ¯ ) + e i | β | 2 x 0 x d t ¯ c 1 e i 1 2 β t ¯ ( | α | t ¯ ) 1 2 W 1 4 i π , 1 4 ( i 1 2 | α | t ¯ 2 ) Wron s 2 ( 1 ) ( t ¯ ) (82)

First, we consider Equation (81). For large x case, that is t ¯ large, Wron s 1 ( 1 ) ( t ¯ ) and Wron s 1 ( 1 ) ( t ¯ ) become

Wron s 1 ( 1 ) ( t ¯ )= W 1 4 + 1 π , 1 4 ( i 1 2 ( | α | t ¯ + | β | | α | ) 2 )i | α | ( i ) ( | α | t ¯ + | β | | α | ) 1 2 (83)

Wron s 1 ( 1 ) ( t ¯ )= W 1 4 + 1 π , 1 4 ( i 1 2 ( | α | t | β | | α | ) 2 )i | α | ( | α | t ¯ | β | | α | ) 1 2 (84)

Then, I 1 ( 1 ) becomes

I 1 ( 1 ) =i e i | β | 2 x x d t ¯ c 1 e i 1 2 α t ¯ 2 e i 1 2 β t ¯ ( | α | t ¯ ) 1 2 W 1 4 i π , 1 4 ( i 1 2 | α | t ¯ 2 ) W 1 4 + 1 π , 1 4 ( i 1 2 ( | α | t ¯ + | β | | α | ) 2 )i | α | ( | α | t ¯ + | β | | α | ) 1 2 + e i | β | 2 x x d t ¯ c 1 e i 1 2 α t ¯ 2 e i 1 2 β t ¯ ( | α | t ¯ ) 1 2 W 1 4 i π , 1 4 ( i 1 2 | α | t ¯ 2 ) W 1 4 + 1 π , 1 4 ( i 1 2 ( | α | t ¯ | β | | α | ) 2 )i | α | ( | α | t ¯ | β | | α | ) 1 2 (85)

Note that we need upper limit of integration range only because we use approximated form of Wronskian for large t ¯ case.

We use asymptotic form of W κ,μ that is defined as [36]

W κ,μ ( z )~ e z 2 z κ [ 1+ n=1 [ μ 2 ( κ 1 2 ) 2 ][ μ 2 ( κn+1 ) 2 ] n! z n ] (86)

Recalling that ( i 1 2 | α | t ¯ 2 ) 1 4 i π = ( 1 2 ) 1 4 i π e 1 2 +i π 8 ( | α | x ) 1 2 2i π , ( i 1 2 ( | α | t ¯ ± | β | | α | ) 2 ) 1 4 + i π = ( 1 2 ) 1 4 + i π e 1 2 +i π 8 ( | α | x± | β | | α | ) 1 2 + 2i π .

I 1 ( 1 ) ~ ( 1 2 ) 1 2 2i π i e i | β | 2 x x d t ¯ e i | β | 2 t ¯ ( | α | t ¯ ) 1 2i π c 1 i | α | 1+ i π t ¯ ( t ¯ + | β | | α | ) 1 2 + 2i π + ( 1 2 ) 1 2 2i π e i | β | 2 x x d t ¯ e i | β | 2 t ¯ ( | α | t ¯ ) 1 2i π c 1 i | α | 1+ i π t ¯ ( t ¯ | β | | α | ) 1 2 + 2i π (87)

Because t ¯ is large and | β | is small, essential integral of Equation (87) is expressed as

x d t ¯ e i | β | 2 t ¯ 1 t ¯ 2+ 4i π and x d t ¯ e i | β | 2 t ¯ 1 t ¯ 2+ 4i π (88)

both integral in Equation (88) become as

x d t ¯ e ±i | β | 2 t ¯ 1 t ¯ 2+ 4i π = e ±i | β | 2 x x 1+ 4i π ( ±i | β | 2 ) x d t ¯ e ±i | β | 2 t ¯ t ¯ 1+ 4i π (89)

Because | β | is small, the first term is sufficient. Then, Equation (87) becomes

I 1 ( 1 ) ~ c ¯ 1 ( 1+i ) 1 x 1+ 4i π (90)

Then, for very small | q | case, Fourier Transform of the first part of particular solution F 3( 1 ) ( 1 ) ( | q | ) becomes

F 3( 1 ) ( 1 ) ( | q | )= 0 dx e i| q |x c ¯ 1 i 1 x 1+ 4i π e i 1 4 | α | x 2 ( | α | x+ | β | | α | ) 1 2 W 1 4 + i π , 1 4 ( i 1 2 ( | α | x+ | β | | α | ) 2 ) + 0 dx e i| q |x c ¯ 1 1 x 1+ 4i π e i 1 4 | α | x 2 ( | α | x | β | | α | ) 1 2 W 1 4 + i π , 1 4 ( i 1 2 ( | α | x | β | | α | ) 2 ) = 0 dx e i| q |x c ¯ 1 ( 1+i ) 1 x 1+ 4i π e i 1 4 | α | x 2 ( | α | x ) 1 2 W 1 4 + i π , 1 4 ( i 1 2 | α | x 2 ) (91)

Note that we obtain the last line of Equation (91) using the condition that | β | is small and x is large.

Changing variable as x= e i π 4 z ¯ | α | and also considering contour integral as shown in ref. [23], apart from constant F 3 ( 1 ) ( | q | ) is expressed as

F 3( 1 ) ( 1 ) ( | q | )= 0 d z ¯ e μ + z ¯ 1 z ¯ 1+ 2i π 0 dt e t t 1 i π ( 1+ 2t z 2 ) 1 2 + i π (92)

where μ + = | q | 2| α | ( 1+i ) .

Note that, for κ= 1 4 + i π case, the term ( 1 2 ( | α | x ) 2 ) 1 4 + i π , come from W 1 4 + i π , 1 4 ( 1 2 ( | α | x ) 2 ) , of which real part is proportional to ( | α | x ) 1 2 cancel out ( | α | x ) 1 2 term.

In Equation (92), we use integral representation of W κ,μ which is defined as [36]

W κ,μ ( z )= e z 2 z κ Γ( μκ+ 1 2 ) 0 dt e t t μκ 1 2 ( 1+ t z ) μ+κ 1 2 (93)

Equation (92) can be evaluated as follows.

F 3( 1 ) ( 1 ) ( | q | )= 0 dt e t t 1 i π 0 d z ¯ e μ + z ¯ 1 z ¯ 1+ 2i π ( z ¯ 2 z ¯ 2 +2t ) 1 2 i π = 0 dt e t t 1 i π 0 d z ¯ e μ + z ¯ z ¯ 4i π ( z ¯ 2 +2t ) 1 2 + i π = 1 π Γ( 1 2 i π ) 0 dt e t t 1 i π G 13 31 ( μ + 2 t 2 | 1 2 + 2i π i π 0 1 2 ) (94)

where G 13 31 denotes Meijer’s G-function defined as [37]

G pq mn ( x| a r b r )= h=1 m ( j=1 m Γ( b j b h ) ) j=1 n Γ( 1+ b h a j ) j=m+1 w Γ( 1+ b h b j ) j=n+1 p Γ( a j b h ) x b h × F p q1 ( 1+ b h a 1 ,,1+ b h a p ;1+ b h b 1 ,,1+ b h b q : ( 1 ) pmn x )

The prime by the product symbol denotes the omission of the product when j = h. The asterisk under the symbol for the function F p q1 denotes the omission of the h-th parameter. This is defined under the condition that either p<q or P=q and | x |<1 .

F p q is called generalized hyper geometric series defined as

F p q ( α 1 ,, α p ; β 1 ,, β k :z )= k=0 ( α 1 ) k ( α p ) k ( β 1 ) k ( β q ) k z k k!

where ( α j ) k = α j ( α j +1 )( α j +k1 ) .

To obtain the last line of Equation (94), we use the following formula [37].

0 dx x 2ν1 ( s 2 + x 2 ) η1 e μx = s 2ν+2η2 2 π Γ( 1η ) G 13 31 ( μ 2 s 2 4 | 1ν 1ην 0 1 2 ) (95)

By following the argument given in ref. [23], Equation (94) shows that the first term of F 3( 1 ) ( 1 ) ( | q | ) is constant or linear of | q | because of μ + = | q | 2| α | ( 1+i ) . We choose the latter case so that the first term of F 3( 1 ) ( 1 ) ( | q | ) is linear of | q | because this satisfies the similarity to our pion wave function. This choice means t term of G-function starts t so that this guarantees that Equation (94) is well defined because t terms come from G-function satisfies the condition that real part of exponent of t terms in integrand becomes larger than −1.

For the other part of F 3 ( 1 ) ( | q | ) denoted as F 3( 2 ) ( 1 ) , we consider Equation (82). Then, for large x case, Wronskia n 2 ( 1 ) ( t ) and Wronskia n 2 ( 1 ) ( t ) become as

Wronskia n 2 ( 1 ) ( t )= W 1 4 1 π , 1 4 ( i 1 2 ( | α | t ¯ + | β | | α | ) 2 )i | α | ( i ) ( | α | t ¯ + | β | | α | ) 1 2 (96)

Wronskia n 2 ( 1 ) ( t )= W 1 4 1 π , 1 4 ( i 1 2 ( | α | t ¯ + | β | | α | ) 2 )i | α | ( | α | t ¯ + | β | | α | ) 1 2 (97)

Then, I 2 ( 1 ) becomes

I 2 ( 1 ) =i e i | β | 2 x x d t ¯ c 1 e i 1 2 β t ¯ ( | α | t ¯ ) 1 2 W 1 4 i π , 1 4 ( i 1 2 | α | t ¯ 2 ) W 1 4 1 π , 1 4 ( i 1 2 ( | α | t ¯ + | β | | α | ) 2 )i | α | ( | α | t ¯ + | β | | α | ) 1 2 + e i | β | 2 x x d t ¯ c 1 e i 1 2 β t ¯ ( | α | t ¯ ) 1 2 W 1 4 i π , 1 4 ( i 1 2 | α | t ¯ 2 ) W 1 4 1 π , 1 4 ( i 1 2 ( | α | t ¯ | β | | α | ) 2 )i | α | ( | α | t ¯ | β | | α | ) 1 2 (98)

Because we are considering in the case of large x case (large t case) and | β | small, we can cancel out Whittaker function. Then, integral of Equation (98) essentially become as

x d t ¯ e i | β | 2 t 1 t ¯ ( t ¯  +  | β | | α | ) 1 2 (99)

x d t ¯ e i | β | 2 t 1 t ¯ ( t ¯ | β | | α | ) 1 2 (100)

To evaluate Equation (99), we employ the following way.

x d t ¯ e i | β | 2 t 1 t ¯ ( t ¯ + | β | | α | ) 1 2 = x d t ¯ e i | β | 2 t ¯ t ¯ ( t ¯ + | β | | α | ) 1 2 t ¯ ( t ¯ + | β | | α | ) = x d t ¯ e i | β | 2 t ¯ 1 | β | | α | ( ( t ¯ + | β | | α | ) 1 2 t ¯ t ¯ ( t ¯ + | β | | α | ) 1 2 ) =2 e i | β | 2 x 1 | β | | α | x ( x+ | β | | α | ) 1 2 i| α | x d t ¯ e i | β | 2 t ¯ t ¯ ( t ¯ + | β | | α | ) 1 2 (101)

2 x d t ¯ e i | β | 2 t ¯ 1 | β | | α | t ¯ ( t ¯ + | β | | α | ) 1 2

Using the condition that t ¯ is large and | β | small, second term of the second line of Equation (101) can be evaluated as

second term=2 x d t ¯ e i | β | 2 t ¯ 1 | β | | α | ( 1 1 2 | β | | α | t ¯ ) =2 | α | | β | e i | β | 2 x xi| α | x d t ¯ e i | β | 2 t ¯ t ¯ x d t ¯ e i | β | 2 t ¯ t ¯ =2 | α | | β | e i | β | 2 x x (102)

Recalling that | β | is small for both Equation (101) and Equation (102), after multiplying e i | β | 2 x term, we obtain

e i | β | 2 x x d t ¯ e i | β | 2 t 1 t ¯ ( t ¯ + | β | | α | ) 1 2 =2 1 | β | | α | x ( x+ | β | | α | ) 1 2 2 1 | β | | α | x =2 1 | β | | α | x( 1+ 1 2 | β | | α | x )2 1 | β | | α | x=1 (103)

Using similar argument for Equation (100), Equation (100) becomes

e i | β | 2 x x d t ¯ e i | β | 2 t 1 t ¯ ( t ¯ | β | | α | ) 1 2 = x d t ¯ e i | β | 2 t ¯ 1 | β | | α | ( t ¯ ( t ¯ | β | | α | ) 1 2 ( t ¯ | β | | α | ) 1 2 t ¯ ) =2 1 | β | | α | x ( x | β | | α | ) 1 2 2 1 | β | | α | x=1 (104)

Because we are considering large x case, we can describe I 2 ( 1 ) in Equation (98) as ( i1 ) .

Then, the second part of Fourier Transform of particular solutions F 3( 2 ) ( 1 ) ( | q | ) for very small | q | case is expressed as

F 3( 2 ) ( 1 ) ( | q | )= 0 dx e i| q |x c ¯ 1 i e i 1 4 | α | x 2 ( | α | x+ | β | | α | ) 1 2 W 1 4 i π , 1 4 ( i 1 2 ( | α | x+ | β | | α | ) 2 ) + 0 dx e i| q |x c ¯ 1 e i 1 4 | α | x 2 ( | α | x | β | | α | ) 1 2 W 1 4 i π , 1 4 ( i 1 2 ( | α | x | β | | α | ) 2 ) = 0 dx e i| q |x c ¯ 1 ( i1 ) e i 1 4 | α | x 2 ( | α | x ) 1 2 W 1 4 i π , 1 4 ( i 1 2 | α | x 2 ) (105)

Changing variable as x= e i 3 4 π z ¯ | α | and considering contour integration as shown in ref. [23] and using integral representation of Whittaker function shown in Equation (93), apart from constants, F 3( 2 ) ( 1 ) ( | q | ) becomes

F 3( 2 ) ( 1 ) ( | q | )= 0 d z ¯ e μ z ¯ z ¯ 1 2i π 0 dt e t t 1 2 + i π ( 1+ 2t z 2 ) 1 i π = 0 dt e t t 1 2 + i π 0 d z ¯ e μ z ¯ z ¯ 1 2i π ( z ¯ 2 z ¯ 2 +2t ) 1+ i π = 0 dt e t t 1 2 + i π 0 d z ¯ e μ z ¯ z ¯ ( z ¯ 2 +2t ) 1 i π = 0 dt e t t 1 2 + 1 π G 13 31 ( μ 2 t 2 | 0 i π 0 1 2 ) (106)

where μ = | q | 2| α | ( 1i ) .

Note that, for κ= 1 4 i π case, the term ( 1 2 ( | α | x ) 2 ) 1 4 i π ( | α | x ) 1 2 2i π come from W 1 4 i π , 1 4 ( 1 2 ( | α | x ) 2 ) multiplies to ( | α | x ) 1 2 , so that this part becomes ( | α | x ) 1 2i π .

Equation (106) shows that the first term of F 3( 2 ) ( 1 ) ( | q | ) is constant or linear of | q | . Again, we choose the case that the first term of F 3( 2 ) ( 1 ) ( | q | ) is linear of | q | .

We can employ same argument in the case of lower half hemisphere group. Then, apart from constants, we obtain

F 3( 1 ) ( 2 ) ( | q | )= 0 dt e t t 1 i π G 13 31 ( μ 2 t 2 | 1 2 2i π i π 0 1 2 ) (107)

F 3( 2 ) ( 2 ) ( | q | )= 0 dt e t t 1 2 + i π G 13 31 ( μ 2 t 2 | 0 i π 0 1 2 ) (108)

Equation (107) and Equation (108) show that, for small | q | case, the first term of Fourier Transform of particular solutions for lower half hemisphere group is also constant or linear of | q | . Thus, we can choose the condition for very small | q | case that the first term of Fourier Transform of particular solutions is linear of | q | .

Because the first term of expansion of Equation (80) becomes 1 | q | , simple way for satisfying the condition that, for very small | q | case, the first term of Fourier Transform of particular solutions is linear of | q | and that, for large | q | case, Fourier Transform of particular solutions asymptotically approaches 1 | q | 3 is multiplying ( 1 e | q | 2| α | u ) 2 to Equation (80).

Then, our form of kaon distribution amplitude | F 3k ( | q | ) | is expressed as

| F 3k ( | q | ) |= 1 | q | 3 ( 1exp( | q |ρ ) ) 2 ( 1exp( | q |ρ )( sin( | q |ρ )+cos( | q |ρ ) ) ) (109)

where ρ= u 2| α | .

Recalling the argument in pion distribution amplitude such that

| q |=( dimension of momentum )×( dimensionless| q | ) =( dimension of momentum )×( dimensionless q ¯ 1 q ¯ )

and (dimension of momentum) moves to u in ρ , and that region of q ¯ is [ 0,1 ] , Equation (107) becomes

| F 3k ( x ) |= ( 1x x ) 3 ( 1exp( x 1x ρ ) ) 2 ×( 1exp( x 1x ρ )( sin( x 1x ρ )+cos( x 1x ρ ) ) ) (110)

Note that we change the notation from q ¯ to x (fraction) [ 0,1 ] in Equation (110) and that actual form of 1 | q | 3 in Equation (109) is 1 ( | q | 2| α | ) 3 that is dimensionless. From now on, x denotes fraction. Because distribution of amplitude is defined by normalized form, we omit the factor of 2| α | in Equation (109).

As we mentioned before, distribution amplitude is defined by normalized form and distribution function is defined by multiplying x to distribution amplitude.

Thus, our kaon valence quark distribution amplitude F 3k ( A ) ( x ) is described as

F 3k ( A ) ( x )= | F 3k ( x ) | 0 1 dx | F 3k ( x ) | (111)

Kaon valence quark distribution functions f k v ( x ) is expressed as

f k v ( x )=x F 3k ( A ) ( x ) (112)

Then, kaon flavor-quark distribution functions becomes

f k ( q ) ( x )= 1 2 x F 3k ( A ) ( x ) (113)

Figure 3 shows comparison of u-quark distribution functions to s-quark distribution function of kaon.

Figure 3. Comparison of u-quark distribution function to s-quark distribution function of kaon: s-quark of kaon (blue curve) and u-quark of kaon (red curve).

Because mass of s-quark is much larger than that of u-quark, scale for evolution of s-quark Q 1 2 is different from that of u-quark Q 2 2 so that we use two valence quark distribution functions (each corresponds different ρ values), i.e., in our case, ρ=0.6 and ρ=1.5 . We choose these ρ values for obtaining appropriate u k / u π ratio shown in later. Because kaon valence quark is u and s ¯ ( K + ) or u ¯ and s ( K ), valence quark distribution functions for both cases are composed by u-quark distribution functions and s-quark distribution functions. Thus, obtained u-quark and s-quark distribution functions are combination of two ρ values valence quark distribution functions. For our case, we determine s-quark distribution functions and u-quark distribution functions as follows.

f k ( s ) ( x )= 1 2 ( d 1 x F 3k( ρ=0.6 ) ( A ) ( x )+ d 2 x F 3k( ρ=1.5 ) ( A ) ( x ) ) (114)

f k ( u ) ( x )= 1 2 ( d 2 x F 3k( ρ=0.6 ) ( A ) ( x )+ d 1 x F 3k( ρ=1.5 ) ( A ) ( x ) ) (115)

where

d 1 = 0 1 dx F 3k ( ρ=0.6 ) ( x ) 0 1 dx F 3k ( ρ=0.6 ) ( x )+ 0 1 dx F 3k ( ρ=1.5 ) ( x )

d 2 = 0 1 dx F 3k ( ρ=1.5 ) ( x ) 0 1 dx F 3k ( ρ=0.6 ) ( x )+ 0 1 dx F 3k ( ρ=1.5 ) ( x )

Figure 4 shows the ratio of s-quark distribution functions of kaon to u-quark distribution functions of kaon. Our results are within 1σ uncertainty region shown in the latest JAM analysis [22].

Figure 4. Ratio of distribution function of s-quark of kaon to that of u-quark of kaon fs(k)/fu(k) (grey curve).

Figure 5 shows our kaon u-quark distribution functions and the ratio of that to pion u-quark distribution functions u ( k ) / u ( π ) . Data denotes experiment results of π-induced Drell-Yan measurements that can be interpreted in terms of K / π structure function ratio which is related to the ratio of u-quark distribution functions of u ( k ) / u ( π ) [18].

Figure 5. Comparison of u-quark distribution functions of kaon to that of pion u-quark distributions of kaon (blue curve), u-quark distributions of pion (red curve), ratio of u-quark distribution functions of kaon to that of pion (grey curve), and data (yellow dots).

For the ratio case, the discrepancy becomes larger in the range of x0.9 . This is because the exponent of ( 1x ) when x approaches 1 is 3 for our distribution functions of kaon, while that of pion is 1.

4. Results and Summary

We obtain pion valence quark distribution functions and pion u-quark distribution functions as shown in Figure 1 and Figure 2, respectively. Exponent of ( 1x ) of asymptotic limit of our pion valence quark and u-quark distribution functions is 1. Our results in Figure 2 follow E615 original analysis data. In this regard, even Dyson Schwinger equation method, adopting inhomogeneous Bethe-Salpeter equation, Bedner et al. [11] show that their results are in excellent agreement with E615 original data over the entire x domain of the data despite imposing ( 1x ) 2 in their representation of pion distribution functions. In addition, very recently, Francis et al. show exponent of ( 1x ) is almost 1 in Lattice QCD calculation [38]. We think that their results are very suggestive. As we mentioned in ref. [23], exponent of asymptotic limit of our pion distribution function, that is 1, corresponds to q 1 in 3 + 1 D by Drell-Yan-West relations in ref. [28]. This is different from current experiment results of charged pion (most likely q 2 ) but is consistent with our results in ref. [24].

In the case of kaon distribution amplitude and functions, our bare kaon distribution amplitude is described as Equation (110) and using normalization and multiplying x, we obtain kaon distribution functions. Our kaon u-quark and s-quark distribution functions are shown in Figure 3. The ratio of s-quark distribution functions to u-quark distribution functions is shown in Figure 4. The ratio of kaon u-quark distribution functions to pion u-quark distribution functions is shown in Figure 5.

Conflicts of Interest

The author declares no conflicts of interest regarding the publication of this paper.

References

[1] Conway, J.S., et al. (1989) Experimental Study of Muon Pairs Produced by 252-GeV Pions on Tungsten. Physical Review D, 39, 92-122.
[2] Aicher, M., Schäfer, A. and Vogelsang, W. (2010) Soft-Gluon Resummation and the Valence Parton Distribution Function of the Pion. Physical Review Letters, 105, Article ID: 252003.[CrossRef] [PubMed]
[3] Frederico, T. and Miller, G.A. (1994) Deep-Inelastic Structure Function of the Pion in the Null-Plane Phenomenology. Physical Review D, 50, 210-216.[CrossRef] [PubMed]
[4] Broniowski, W., Arriola, E.R. and Golec-Biernat, K. (2008) Generalized Parton Distributions of the Pion in Chiral Quark Models and Their QCD Evolution. Physical Review D, 77, Article ID: 034023.[CrossRef]
[5] Melnitchouk, W. (2003) Quark-Hadron Duality in Electron-Pion Scattering. The European Physical Journal A, 17, 223-234.[CrossRef]
[6] Pasquini, B., Rodini, S. and Venturini, S. (2023) Valence Quark, Sea, and Gluon Content of the Pion from the Parton Distribution Functions and the Electromagnetic Form Factor. Physical Review D, 107, Article ID: 114023.[CrossRef]
[7] Xie, G., Li, M., Han, C., Wang, R. and Chen, X. (2021) Simulation of Neutron-Tagged Deep Inelastic Scattering at EICC. Chinese Physics C, 45, Article ID: 053002.[CrossRef]
[8] Nguyen, T., Bashir, A., Roberts, C.D. and Tandy, P.C. (2011) Pion and Kaon Valence-Quark Parton Distribution Functions. Physical Review C, 83, Article ID: 062201.[CrossRef]
[9] Chang, L. and Thomas, A.W. (2015) Pion Valence-Quark Parton Distribution Function. Physics Letters B, 749, 547-550.[CrossRef]
[10] Shi, C., Mezrag, C. and Zong, H. (2018) Pion and Kaon Valence Quark Distribution Functions from Dyson-Schwinger Equations. Physical Review D, 98, Article ID: 054029.[CrossRef]
[11] Bednar, K.D., Cloët, I.C. and Tandy, P.C. (2020) Distinguishing Quarks and Gluons in Pion and Kaon Parton Distribution Functions. Physical Review Letters, 124, Article ID: 042002.[CrossRef] [PubMed]
[12] Ding, M., Raya, K., Binosi, D., Chang, L., Roberts, C.D. and Schmidt, S.M. (2020) Drawing Insights from Pion Parton Distributions. Chinese Physics C, 44, Article ID: 031002.[CrossRef]
[13] Gao, X., Hanlon, A.D., Karthik, N., Mukherjee, S., Petreczky, P., Scior, P., et al. (2022) Continuum-Extrapolated NNLO Valence PDF of the Pion at the Physical Point. Physical Review D, 106, Article ID: 114510.[CrossRef]
[14] Lu, Y., Xu, Y., Raya, K., Roberts, C.D. and Rodríguez-Quintero, J. (2024) Pion Distribution Functions from Low-Order Mellin Moments. Physics Letters B, 850, Article ID: 138534.[CrossRef]
[15] Barry, P.C., Sato, N., Melnitchouk, W. and Ji, C. (2018) First Monte Carlo Global QCD Analysis of Pion Parton Distributions. Physical Review Letters, 121, Article ID: 152001.[CrossRef] [PubMed]
[16] Novikov, I., Abdolmaleki, H., Britzger, D., Cooper-Sarkar, A., Giuli, F., Glazov, A., et al. (2020) Parton Distribution Functions of the Charged Pion within the Xfitter Framework. Physical Review D, 102, Article ID: 014040.[CrossRef]
[17] Alexandrou, C., Bacchio, S., Constantinou, M., Delmar, J., Finkenrath, J., Kostrzewa, B., et al. (2025) Quark and Gluon Momentum Fractions in the Pion and in the Kaon. Physical Review Letters, 134, Article ID: 131902.[CrossRef] [PubMed]
[18] Badier, J., Boucrot, J., Bourotte, J., Burgun, G., Callot, O., Charpentier, Ph., et al. (1980) Measurement of the K/π Structure Function Ratio Using Drell-Yan Process. Physics Letters B, 93, 354-356.[CrossRef]
[19] Badier, J., Boucrot, J., Bourotte, J., Burgun, G., Callot, O., Charpentier, Ph., et al. (1983) Experimental J/Ψ Hadronic Production from 150 to 280 GeV/c. Zeitschrift für Physik C Particles and Fields, 20, 101-116.[CrossRef]
[20] Xu, Z., Binosi, D., Chen, C., Raya, K., Roberts, C.D. and Rodríguez-Quintero, J. (2025) Kaon Distribution Functions from Empirical Information. Physics Letters B, 865, Article ID: 139451.[CrossRef]
[21] Bourrely, C., Buccella, F., Chang, W. and Peng, J. (2024) Extraction of Kaon Partonic Distribution Functions from Drell-Yan and J/ψ Production Data. Physics Letters B, 848, Article ID: 138395.[CrossRef]
[22] Barry, P.C., Ji, C.R., Melnitchouk, W., Sato, N. and Steffens, F. (2025) First Simultaneous Global QCD Analysis of Kaon and Pion Parton Distributions with Lattice QCD Constraints. arXiv: 2510.119v1.
[23] Kurai, T. (2025) New Approach to Pion Distribution Amplitude. Journal of Modern Physics, 16, 886-910.[CrossRef]
[24] Kurai, T. (2021) Light Meson Mass Spectra with Massive Quarks. Journal of Modern Physics, 12, 1545-1572.[CrossRef]
[25] Lu, Y., Chang, L., Raya, K., Roberts, C.D. and Rodríguez-Quintero, J. (2022) Proton and Pion Distribution Functions in Counterpoint. Physics Letters B, 830, Article ID: 137130.[CrossRef]
[26] Miyama, M. and Kumano, S. (1996) Numerical Solutions of Q2 Evolution Equations in a Brute-Force Method. Computer Physics Communications, 94, 185-215.[CrossRef]
[27] Wu, Q., Han, C., Qing, D., Kou, W., Chen, X., Wang, F., et al. (2023) Pion Parton Distribution Functions with the Nonrelativistic Constituent Quark Model. Nuclear Physics B, 994, Article ID: 116321.[CrossRef]
[28] Arrington, J., Ayerbe Gayoso, C., Barry, P.C., Berdnikov, V., Binosi, D., Chang, L., et al. (2021) Revealing the Structure of Light Pseudoscalar Mesons at the Electron-Ion Collider. Journal of Physics G: Nuclear and Particle Physics, 48, Article ID: 075106.[CrossRef]
[29] Kurai, T. (2014) The Meson as a Bound System in 2D Quantum Chromodynamics. Progress of Theoretical and Experimental Physics, 2014, 53B01.[CrossRef]
[30] Suura, H. (1978) Derivation of a Quark-Confinement Equation in the Hamiltonian Formalism of Gauge Field Theories. Physical Review D, 17, 469-482.[CrossRef]
[31] Kurai, T. (2018) Light Meson Mass Spectra and Pion Electromagnetic Form Factor as a Bound System in 3 + 1 Dimensional QCD. Results in Physics, 10, 865-881.[CrossRef]
[32] Weinberg, S. (2000) The Quantum Theory of Field III. Cambridge University Press.
[33] Casher, A., Kogut, J. and Susskind, L. (1974) Vacuum Polarization and the Absence of Free Quarks. Physical Review D, 10, 732-745.[CrossRef]
[34] Gakhov, F.D. (1966) Riemann Boundary Value Problem. In: Gakhov, F.D., Ed., Boundary Value Problems, Elsevier, 85-142.[CrossRef]
[35] Ince, E.L. (1956) Ordinary Differential Equations. Dover Publication.
[36] Moriguchi, S., Udagawa, K. and Hitotsumatsu, S. (1975) Mathematics Formula III. Iwanami.
[37] Gradshteyn, I.S. and Ryzhik, M. (1980) Table of Integrals, Series, and Products (Corrected and Enlarged Edition). Academic Press.
[38] Francis, A., Fritzsch, P., Karur, R., Kim, J., Pederiva, G., Pefcou, D.A., et al. (2025) Moments of Parton Distributions Functions of the Pion from Lattice QCD Using Gradient Flow. arXiv: 2510.26738.

Copyright © 2026 by authors and Scientific Research Publishing Inc.

Creative Commons License

This work and the related PDF file are licensed under a Creative Commons Attribution 4.0 International License.