Extracting the longitudinal structure function FL(x,Q2) at small x from a Froissart-bounded parametrization of F2(x,Q2)

We present a method to extract, in the leading and next-to-leading order approximations, the longitudinal deep-inelastic scattering structure function FL(x,Q2) from the experimental data by relying on a Froissart-bounded parametrization of the transversal structure function F2(x,Q2) and, partially, on the Dokshitzer-Gribov-Lipatov-Altarelli-Parisi equations. Particular attention is paid on kinematics of low and ultra low values of the Bjorken variable x. Analytical expressions for FL(x,Q2) in terms of the effective parameters of the parametrization of F2(x,Q2) are presented explicitly. We argue that the obtained structure functions FL(x,Q2) within both, the leading and next-to-leading order approximations, manifestly obey the Froissart boundary conditions. Numerical calculations and comparison with available data from ZEUS and H1-Collaborations at HERA demonstrate that the suggested method provides reliable structure functions FL(x,Q2) at low x in a wide range of the momentum transfer (1 GeV2


I. INTRODUCTION
At small values of the Bjorken variable x, the nonperturbative effects in the deep inelastic structure functions (SF) were expected to play a decisive role in describing the corresponding cross sections. However, it has been observed, cf. Ref. [1], that even in the region of low momentum transfer Q 2 ∼ 1 GeV 2 , where traditionally the soft processes were considered to govern the cross sections, the perturbative QCD (pQCD) methods could be still adequate in description of high energy processes, in particular, at relatively low values of x: 10 −5 ≤ x ≤ 10 −2 . It should also be noted that, at extremely low x, x → 0, the pQCD evolution leads, nonetheless, to a rather singular behaviour of the parton distribution functions (PDF) (see e.g. Ref. [2] and references therein quoted), which is in a strong disagreement with the Froissart boundary conditions [3]. In Refs. [4][5][6] M. M. Block et al. have suggested a new parametrization of the SF F 2 (x, Q 2 ) which describes fairly well the available experimental data on the reduced cross sections and, at asymptotically low x, provides a behavior of the hadron-hadron cross sections ∼ ln 2 s at large s, where s is the Mandelstam variable denoting the square of the total invariant energy of the process, in a full accordance with the Froissart predictions [3]. The most recent parametrization suggested in Ref. [6] by M. M. Block, L. Durand and P. Ha, in what follows referred to as the BDH parametrization, is also pertinent in investigations of leptonhadron processes at ultra-high energies , e.g. the scattering of cosmic neutrinos from hadrons [5][6][7][8][9]. Note that, in case of neutrino scattering other SF's, such as the pure valence, F 3 (x, Q 2 ), and longitudinal, F L (x, Q 2 ), are relevant to describe the process. While at low values of x the valence structure function F 3 (x, Q 2 ) vanishes, the longitudinal F L (x, Q 2 ) remains finite and can even be predominant in the cross section. Thus, a theoretical analysis of the longitudinal SF F L (x, Q 2 ) at low x, in context of fulfilment of the Froissart prescriptions, is of a great importance in treatments of ultra-high energy processes as well.
Hitherto, most theoretical analyses [8,10] of neutrino processes have been performed in the leading order (LO) approximation, within which the Callan-Gross relation is assumed to be satisfied exactly, i.e. the longitudinal structure function F L = 0. Beyond the LO the effects from F L can be sizable, hence it can not be longer neglected, cf. Ref. [8,11].
In the present paper we present a method of extraction of the longitudinal SF, F L (x, Q 2 ) in the kinematical region of low values of the Bjorken variable x from the known structure function F BDH 2 (x, Q 2 ) and known derivative dF BDH 2 /d ln(Q 2 ) by relying, with some extent, on the Dokshitzer-Gribov-Lipatov-Altarelli-Parisi (DGLAP) Q 2 -evolution equations [12]. In our calculations we use the most recent version of the BDH-parametrization reported in Ref. [6]. In fact, the presented approach is a further development of the methods previously suggested in Refs. [13,14] to extract some general characteristics of the gluon density and longitudinal SF at low x from the the experimentally known SF F 2 (x, Q 2 ) and logarithmic derivative dF 2 /d ln Q 2 .
The extraction procedure has been inspired by the Altarelli-Martinelli formula [15] suggested to determine the gluon density from F L (x, Q 2 ), and improved in Ref. [16].
In our case, the SF F 2 (x, Q 2 ) is considered experimentally known being defined by the F BDH 2 (x, Q 2 ) parametrization [6], i.e. the x and Q 2 dependencies of the transverse SF and the corresponding logarithmic derivative are supposed to be known. As a first step of the analysis, the method has been applied to extract F L (x, Q 2 ) in the LO. Results of such a procedure have been briefly reported in Ref. [17], where it has been demonstrated that the extracted structure function F BDH L (x, Q 2 ) at moderate and low values of x is in a reasonable good agreement with the available experimental data [19]. However, for ultra low values of x the agreement becomes less satisfactorily and even rather poor in the limit x → 0. This serves as a clear indication that the LO analysis is not sufficient in the region x → 0 and the next-to-leading order (NLO) corrections become significant and are to be implemented in to the extraction procedure. Similar investigations of the longitudinal SF have been performed in Ref. [18].
In this paper, we present in some details the LO analysis [17], and provide further development of the method extending it beyond the LO approximation by considering and resumming the NLO corrections.
It is worth emphasizing that the NLO approximation for F L (x, Q 2 ), i.e. calculations up to α 2 s -corrections, corresponds to the next-next-to-leading order (NNLO) approximation for F 2 (x, Q 2 ) which, in the LO, is ∝ α 0 s . Hence, in our approach it becomes possible to perform NLO and NNLO analyses of the ultra-high energy ( √ s ∼ 1 TeV) neutrino cross-sections similar to existing NLO [20] and NNLO [21] investigations based on pQCD. Such analyses are rather important in view of the recently appeared possibility of a direct comparison with emerging data from the IceCube Collaboration [22] (cf. also Ref. [23]) and anticipated data from the IceCube-Gen2 [24], whose performance is much better and which is awaited to provide substantially more precise measurements of the neutrino-nucleon cross-section.
Our paper is organized as follows: In Sec. II we present the basic formulae of the approach.

II. BASIC FORMULAE
In view of at low values of x the non-singlet quark distributions become negligibly small in comparison with the singlet distributions, in the present analysis they are disregarded. Then, the transverse F 2 (x, Q 2 ) and longitudinal F L (x, Q 2 ) structure functions are expressed solely via the singlet quark and gluon densities xf a (x, Q 2 ) (hereafter a = s, g and k = 2, L) as where e is the average charge squared, f with f as the number of considered flavors, q 2 = −Q 2 and x = Q 2 /2pq (p being the momentum of the nucleon) denote the momentum transfer and the Bjorken scaling variable, respectively. The quantities B k,a (x) are the known Wilson coefficient functions. In Eq. (1) and throughout the rest of the paper, the symbol ⊗ is used for a shorthand notation of the convolution formula, i.e.
According to the DGLAP Q 2 -evolution equations [12] the leading twist quark, xf s (x, Q 2 ), and gluon, xf g (x, Q 2 ), distributions obey the following system of integro-differential equations where P ab (x) (a, b = s, g) are the corresponding splitting functions.
Within the pQCD, and up to the NLO corrections, the coefficient functions B k,a (x) and the splitting functions P ab (x) read as where a s (Q 2 ) = α s (Q 2 )/4π is the QCD running coupling, which, for convenience, includes in to its definition an additional factor of 4π in comparison with the standard notation. In the above equations and hereafter the superscripts (0, 1) mark the corresponding order of the perturbation where, for brevity, the following notations have been employed with β 0 and β 1 as the first two coefficients of the QCD β-function In Eq. Few remarks are in order here. As known [25,26], equation (7) in its actual form leads to a too singular behaviour of the gluon distribution at small x, violating the Froissart boundary restrictions. One can go beyond the perturbative theory and try to cure the problem by adding in the r.h.s. of Eq. (7) terms proportional to (xf g ) 2 which make the distribution (hence, the corresponding cross-sections [27]) less singular at the origin and can, in principle, reconcile it with the Froissart requirements. A detailed inspection of Eq. (7) in context of implementation of additional modifications to fulfill the Froissart conditions is beyond the scope of the present paper and in what follows we omit it in our analysis. However, the gluon distribution originating from the omitted Eq. (7) and entering in to the remaining equations (8) and (9) is supposed to have the correct asymptotic behavior, i.e. to be of the same LO form as the BDH parametrization of the F BDH 2 (x, Q 2 ). This conjectures has been confirmed in previous analysis [28] where the early parametrization of F 2 (x, Q 2 ) [5] has been employed to determine the gluon density within the LO. Then, Eq. (8) with the known F BDH 2 (x, Q 2 ), can be considered as the definition of the gluon density xf BDH g (x, Q 2 ) in the whole kinematical interval, cf. Ref. [28].
Consequently, in the system (8)-(9) of two equations with two unknown distributions one can eliminate the gluon part and solve the remaining equation with respect to the longitudinal F L (x, Q 2 ) and express it via the known paramterization of F 2 (x, Q 2 ). With these statements, now we are in a position to solve Eqs. (8) and (9) and to extract the desired longitudinal SF.
Notice that, albeit at the first glance the above equations are relatively simple, direct solving of (8)-(9) actually turns out to be a rather complicate and cumbersome procedure. One can substantially simplify the calculations by considering Eqs. (8)- (9) in the space of Mellin momenta, and taking advantage of the fact the convolution form f 1 (x) ⊗ f 2 (x) in x space becomes merely a product of individual Mellin transforms of the corresponding functions in the space of Mellin momenta. Consequently, all our further calculations we perform in Mellin space.

A. Mellin transforms
The Mellin transform of the PDF's entering in to Eqs. (8)-(9) are defined as where, as before, a, b = s, g and k = 2, L. Then, after some algebra, the Mellin transforms of where the anomalous dimensions γ L,a (n) (i = 0, 1) arẽ The explicit expressions for the anomalous dimensions γ L,a (n) can be found in Ref. [29], whereas the NLO part, B L,N S (n) and B (1) L,a (n), has been reported, for the first time, in Refs. [30,31]. Unfortunately, there are several misprints in the above mentioned references. Erratum for Refs. [30,31] can be found in, e.g. Ref. [33]. Yet, in Ref. [29] a factor two in the LO coefficients B (0) L,a (n) is missed. Observe that as mentioned above, the Mellin transform significantly simplifies our calculations for, the complicated integrodifferential system of equations (8)-(9) in x-space, is translated in to a relatively simple system (15)-(16) of pure algebraic equations. Now, we solve the system (15)-(16) w.r.t the longitudinal Mellin momentum M L (n, Q 2 ) and express it through known momentum M 2 (x, Q 2 ) and the

B. Anomalous dimensions and coefficient functions
Here below we present the explicit expressions for the LO ingredients only. The corresponding expressions for the NLO corrections are rather cumbersome and, as already mentioned, can be found in Refs. [29][30][31]33]. Explicitly, the LO anomalous dimensions γ L,a (n), are as follows: In Eqs. (20)-(21) we introduce, and shall widely use throughout the rest of the paper, the notion of the so-called nested sums, defined as Note that in previous calculations [29][30][31] Coming back to the Mellin transforms (15), (16) and (19), we recall that we are interested in investigation of the PDF's in the region of low x. In the Mellin space it corresponds to small momenta, and at extremely low x, it suffices to restrict the analysis with the first momentum and to study the solutions of (15), (16) and (19) for n = 1 + ω at ω → 0. The nested sums S m (n) (here n is not mandatorily an integer) at positive indices m can be expressed via the familiar Riemann Ψ-functions Ψ(n) as where ζ m are Euler constants and the last series in the r.h.s. of (24) where the above series are well defined for any n including noninteger values.
In what follows we are interested in quantities which contribute to NLO at low x, i.e. in the initial series S −2 (n), S −3 (n) and S −2,1 (n), anomalous dimensions and the coefficient functions, at n = 1 + ω. In the vicinity of ω = 0 we have

III. BDH-LIKE RESULTS
We reiterate that, our analysis is based on Eqs. (8)- (9) or, equivalently, on Eqs. (15), (16) and (19), where the structure function F 2 (x, Q 2 ) is supposed to be known and determined from the existing experimental data. In the present paper we employ the BDH parametrization [6], obtained from a combined fit of the H1 and ZEUS collaborations data [19] in a range of the kinematical variables x and Q 2 , x < 0.01 and 0.15 GeV 2 < Q 2 < 3000 GeV 2 . The explicit expression for the BDH parametrization reads as where the dependence on the effective parameters is encoded in D(Q 2 ) and A m (Q 2 ) with the logarithmic terms L as The performed fit of the experimental data [19] provided the following values of the effective parameters [6]: and a 00 = 0.255 ± 0.016, a 01 · 10 1 = 1.475 ± 0.3025, a 10 · 10 4 = 8.205 ± 4.620, a 11 · 10 2 = −5.148 ± 0.819, a 12 · 10 3 = −4.725 ± 1.010, a 20 · 10 3 = 2.217 ± 0.142, a 21 · 10 2 = 1.244 ± 0.0860, a 22 · 10 4 = 5.958 ± 2.320 . (31) Notice that, the BDH parametrization (27) is written in the (x − Q 2 ) space, whereas the equation (19) is in space of Mellin momenta. Hence, before proceeding with consideration of Eq. (19), we transform the BDH parametrization in to Mellin space. Then, in the transformed parametrization M 2 (n, Q 2 ) we consider the first momenta n = 1 + ω and take the limit ω → 0, which corresponds to low x, and obtain where the quantity P m (ω, ν, L 1 ) stands for the approximate expression of the integral The explicit expression for the integral (33) is relegated in Appendix A. It should be noted that P k (ω, ν, L 1 ), besides the finite part at ω → 0, contains also negative powers of ω, i.e.
is a singular function at ω = 0. Namely these singularities are of our interests in further procedure of solving Eqs. (8)- (19). The strategy is as follow: we disregard the finite part of where P sing Note that P sing which are widely used subsequently.
Noticing that we can write Collecting all results together, we have The next step is to consider Eq. (44) at n = 1 + ω and to perform the series expansion of the anomalous dimensions and the coefficient functions about ω = 0. Restricting the expansion up to terms ∝ ω 2 we obtain Observe that, in the above equation (44) both, the momentum M 2 (x, Q 2 ) and the derivative dM 2 (n, Q 2 )/dlnQ 2 are of a similar form, cf. Eqs. (34) and (41), so that the r.h.s. of (44) with the expansions (46) can be represented as This substantially simplifies the calculations since, in such a case, we can apply the recurrent relations (36) to find the desired coefficients. For instance, for any (known) series of the form one easily obtains the result where It is seen that it is sufficient to determine the coefficients T 0 , T 1 and T 2 to find the solution explicitly read as It is worth mentioning once more that, the above results for C i (i = 0, 1, 2) have been obtained in a straight way avoiding direct comparisons of singularities in (44), as previously reported in [17,28]. Another important result is that the adopted BDH parametrization for the F 2 (x, Q 2 ) led directly to the same BDH-like form of the longitudinal structure functions F BDH L,LO (x, Q 2 ) which, consequently, also obey the Froissart boundary condition where ν L = ν and, according to the quark counting rules [35].

V. NLO ANALYSIS
In this section we discuss the longitudinal structure function within the NLO approximation.
As before, particular attention is paid to the region of asymptotically small x → 0 in context with the Froissart boundary.
Multiplying both sides of Eq. (19) by the factor 1 + a s (Q 2 ) δ where Using equations (44) and (53) where the LO momentum M L,LO (n, Q 2 ) has already been calculated in the previous section, cf.

Eq. (50).
Prior to proceed with the inverse Mellin transforms, it is convenient to extract the singular structure of the NLO coefficients δ L,g (n) andB L,s (n). We have L,s (n) =B with (ω → 0) Then, the equation (56) can be rewritten in the following form Now the inverse Mellin transforms of the last equations can be easily performed (see also Appendix B). The result is where With the considered accuracy the obtained equation (61) can be rewritten as Eventually, the final expression for the longitudinal SFF BDH L (x, Q 2 ) reads as This is our final expression for the longitudinal SF F BDH L (x, Q 2 ) within the NLO approximation for low values of x.

VI. RESULTS
With the explicit form of the basic expressions laid above, we can proceed to extract the longitudinal structure function F L (x, Q 2 ) from data mediated by the BDH-parametrization of F BDH 2 (x, Q 2 ). In our calculations we employ the standard representation for QCD couplings in the LO and NLO (within the MS-scheme) approximations a s (Q 2 ) = 1 The QCD parameter Λ has been extracted from the running coupling α s normalized at the Z-boson mass, α s (M 2 Z ), using the b-and c-quarks thresholds according to Ref. [37]. Applying this procedure to ZEUS data, with α s (M 2 Z ) = 0.1166 [38], we obtain the following results for Λ, cf. Ref. [36] LO  negative and result in a better agreement with data. However, at extremely low momenta, Q 2 < 1.5 GeV 2 , the extracted SF within NLO is still above the experimental data. It should also be mentioned that our calculations are consistent with other theoretical results, obtained, e.g. in the framework of perturbation theory [40,41] and/or in the Ref. [42] in the so-called k t -factorization approach [43], both of which incorporate the Balitsky-Fadin-Kuraev-Lipatov (BFKL) resummation [44] at low x (for a review of low-x phenomenology see, e.g. cf. Ref. [45]).
Recall that inclusion of the BFKL ressummation in study of PDF's lead to an improvement of the description of data at small x and nowadays appears as an integral part in a bulk of approaches. The NNPDF [46] Collaboration and the xFitter HERAPDF team [19,47], whose approaches are based on the DGLAP equations, recently included the BFKL resummation in to their analysis of the combined H1&ZEUS inclusive cross-section [48] achieving, in such a way, a much better description of data [40,41]. Analogous studies have been performed in Refs. [40,[49][50][51]. This is in some contrast with results of the standard PDF sets [46,47,52] without the BFKL resummations.
A particular interests present the ratio of the longitudinal to transversal cross sections, defined as Recently, the H1-Collaboration has reported the ratio R L (x, Q 2 ) measured in several kinematical bins of averaged Q 2 and x, cf. Table 6 of Ref. [39]. Within such kinematics, the invariant mass W changes from W ∼ 230 GeV to W ∼ 184 GeV with increase of Q 2 and x in the selected bins. In Figure 2 we present the ratio (67), calculated with the extracted F BDH L (x, Q 2 ) and  Figure 3 where the x-evolution of F L (x, Q 2 ) is clearly exhibited. It is seen that, for all values of the presented Q 2 , the extracted SF within the NLO approximation is in a much better agreement with data. This persuades us that the obtained SF in the NLO approximation can be pertinent in future analysis of the ultra-high energy neutrino data.
In Figure 4 we present the ratio (67) calculated for the same kinematics as in Figure 3. As in previous calculations, the NLO results are in a better agreement with data. It is also seen from Figure 4 that the NLO ratio R L (x, Q 2 ) exhibits a tendency to be almost independent on x in each bin of Q 2 , decreasing, however, with Q 2 increase, as it should be. We also mention that, as in the case of Q 2 -dependence, our extracted longitudinal SF as a function of x is in a reasonable good agreement with other theoretical predictions, see e.g. [40,41] for pQCD results and/or Ref. [53] for results obtained within the k t -factorization approach. Likewise our results are in a good agreement with the previous investigations reported in Ref. [14], were a similar analysis has been performed within the framework of pQCD, with the experimental data for the transverse SF F 2 (x, Q 2 ) and the logarithmic derivative dF 2 /d ln(Q 2 ).

VII. CONCLUSIONS AND OUTLOOK
In this paper, we present a further development of the method of extraction of the longitudinal DIS structure function F L (x, Q 2 ) suggested in Refs. [13,14,17].  (27) and is presented in Eq. (64). Some comments are in order here. Observe that, the final expression (64) contains the dominant logarithmic terms ∼ ln(1/x) in both the numerator and denominator parts. In principle, due to smallness of the running coupling α s (x, Q 2 ), the denominator part with L C can be rewritten in the numerator as an alternate series leading to a behaviour similar to one known within the pQCD, where the NLO corrections, in the considered kinematical region, are negative and large, while NNLO contributions are positive and also large (see, e.g. Ref. [54]). A resummation of these contribution would allow to avoid such an alternate behaviour. In our case, this is achieved by keeping the logarithmic term in the denominator without expanding it in to series relative to α s (x, Q 2 ). With some extent, our representation of the basic NLO corrections, Eq. (64), can be considered as an effective resummation of the most important, at low x, logarithmic terms in each order of the perturbation theory. Another observation is that by keeping the logarithmic L C terms in the denominator we also manifestly demonstrate that the SF F L (x, Q 2 ) obeys the Froissart conditions. As shown in Appendix B, the logarithmic part in Eq. (64) is a direct consequence of the inverse Mellin transforms of terms ∼ 1/ω.
We have applied the developed method to extract the longitudinal SF within the kinematical conditions corresponding to that available at the HERA collider. It has been found that, at relatively large Q 2 > 10 GeV 2 both, LO and NLO results reproduce fairly well the experimental data. At smaller Q 2 the LO approximation fails to describe data, being systematically well above. Accounting for NLO corrections, which at low x turn out to be negative, substantially improve the description of the SF and the ratio of the longitudinal to transversal cross sections.
However, at extremely low momentum transfer Q 2 1 GeV 2 , the extracted SF still exceeds data.
We have performed an analysis of the x-evolution of the extracted SF. It has been demonstrated that the x-dependence of F L (x, Q 2 ) also reproduces the behavior of the experimental data at low x. Both, the extracted SF and the ratio R L (x, Q 2 ) as functions of x, are in a fairly good agreement with data, herewith the NLO results are in a much better agreement not only with data, but also with other existing theoretical investigations based on perturbative QCD, improved by the BFKL resummation (see Ref. [40,41] and references therein quoted), as well as with the results [42] obtained in the framework of the k t -factorization method [43], also based on the BFKL approach [44]. Calculations of the longitudinal SF based on the traditional pQCD without such improvements turn out to be rather unstable, due to the fact that the subsequent perturbative corrections can be even larger than the previous ones [54,55]. Incorporation of corrections inspired by the BFKL resummation lead to a substantially more stable results for F L (x, Q 2 ) [55] (see also similar investigations in Refs. [40,[49][50][51]), aimed to achieve a good description of the combined H1&ZEUS inclusive cross-sections [48].
Apart from the study of the Froissart boundary restrictions, the knowledge of the F L (x, Q 2 ) at low x is of a great interest in connection with the theoretical treatments of the ultrahigh energy processes with cosmic neutrinos. As already mentioned in Introduction, the NLO approximation for F L (x, Q 2 ), i.e. calculations up to α 2 s -corrections, corresponds to next-nextto-leading order (NNLO) for F 2 (x, Q 2 ) which, in the LO, is ∝ α 0 s . Consequently, with the NLO results (64), in our approach it becomes possible to perform NLO and NNLO analyses of the ultra-high energy ( √ s ∼ 1 TeV) neutrino cross-sections similar to NLO [20] and NNLO [21] investigations based on pQCD. Such calculations are of a great importance in view of expected reliable cross sections from existing and forthcoming data at the IceCube [23] and from the substantially improved IceCube-Gen2 [24] Collaborations. Therefore, a direct comparison of the theoretical predictions with experimental data becomes feasible.
In the kinematical region where the gluon contributions are sizable the (large) corrections within the traditional pQCD can be strongly reduced by a proper change of the factorization and renormalization scales [56]. An analysis of the precise H1&ZEUS combined data [19] obtained within the kinematics near the limit of applicability of pQCD has shown [57] that an employ of effective scales with large parameters provides much smaller high-order perturbative corrections. In such a case the strong couplings decrease as well and, as a rule, calculations with effective scales lead to a better agreement with data (cf. investigations of the NLO corrections in context with high-energy asymptotics of virtual photon-photon collision [58] and studies of F L (x, Q 2 ) in the framework of the k t -fragmentation approach [42]). This encourages us to continue our low-x analysis of the SF's by implementing special change of the factorization and renormalization scales [56]. This is the subject of our further investigations and results will be presented elsewhere.
Furthermore, we plan to improve our approach by accommodating the method to extract, in the LO and NLO approximations, also the gluon densities. We shall note that, the gluon distribution is by far less known, experimentally and theoretically. Even the shape of the gluon density are often taken quite different in different PDF sets [52,59], although considered within the same prerequisites of pQCD. However, the range of variation of gluon density strongly decreases when BFKL resummation is included in the analyses at low x (see the most recent publication [21] and discussion therein).
An extraction of the gluon distributions from experimental data, performed within an approach similar to the one suggested in the present paper, can provide valuable additional information on the problem. Such an analysis can be accomplished by employing the charm, F cc 2 (x, Q 2 ), and beauty, F bb 2 (x, Q 2 ), components of the SF F 2 (x, Q 2 ), which are directly related to the gluon density in the photo-gluon fusion reactions (see Ref. [60] and discussion therein). The extracted SF's can be compared with the recently obtained combined data from the H1&ZEUS-Collaboration [61] for the F cc 2 (x, Q 2 ) and F bb 2 (x, Q 2 ) and with the theoretical predictions [21] based on pQCD with BFKL corrections included, and also with the results [62] obtained in the framework of the k t -fragmentation. Furthermore, the charmed parts of transverse F cc 2 (x, Q 2 ) and longitudinal, F cc L (x, Q 2 ), structure functions being calculated within our approach, can be used to predict the charmed part of the neutrino-nucleon cross-sections at ultra-high energy and to compare with other calculations [21] based on pQCD with BFKL corrections. Yet, the BDH-gluon density itself can serve as a useful tool for estimations of the cross sections with cosmic rays, cf. Ref. [63] (for a most recent review on the subject, see Ref. [64] and references therein quoted). Investigations in this direction are in progress.
In summary, we present a theoretical method to extract, from the experimental data, the longitudinal DIS structure function F L (x, Q 2 ) at low x within the leading and next-to-leading order approximations. Explicit, analytical expressions for the structure function in both, LO and NLO, approximations are obtained in terms of the effective parameters of the Froissartbound parametrization of F 2 (x, Q 2 ) and results of numerical calculations as well as comparisons with available experimental are presented.