*
* ..File: xpns2n.f   UPDATE 7/2000 
*
*
* ..Parametrization of the 3-loop MS(bar) splitting functions P_NS^(2) 
*    for the evolution of unpolarized non-singlet partons, mu_r = mu_f. 
*    The expansion parameter is alpha_s/(4 pi).
*
* ..The two sets spanning the error estimate are called via  IMODN = 1
*    and  IMODN = 2.  Any other value of IMODN invokes their average.
*
* ..The distributions (in the mathematical sense) are given as in eq.
*    (B.27) of Floratos, Kounnas, Lacaze: Nucl. Phys. B192 (1981) 417.
*    The name-endings A, B, and C of the functions below correspond to 
*    the kernel superscripts [2], [3], and [1] in that equation.
*
* ..The P^+ results are based on the lowest six even-integer moments
*    calculated by Larin et al. and  Retey and Vermaseren. The P^- and
*    P^S = P^V - P^-) approximations use the seven lowest odd-integer 
*    moments computed also by Retey and Vermaseren. Also used for P^+
*    and P^- are the leading small-x terms (ln^4 x) as obtained by 
*    Blumlein and Vogt. The exact N_f^2 term has been derived by Gracey.
*
* ..Reference: W.L. van Neerven and A. Vogt, 
*              hep-ph/9907472 and hep-ph/0007362
*    It is appropriate to cite also the above-mentioned sources of
*    information, i.e, ref. [2-4] and [6,7] of hep-ph/0007362.
*
*
* =====================================================================
*
*
* ..This is the regular piece of P^(2)+. 
*
       FUNCTION P2NSPA (Y, NF, IMODN)
*
       IMPLICIT REAL*8 (A-Z)
       INTEGER IMODN, NF
       DL  = LOG (Y)
       DL1 = LOG (1.-Y)
*
       G2NSPA1 = - 1047.590 * DL1 + 843.884 * Y**2  + 98.65 * (1.-Y) 
     2           + 33.71 * DL**2  - 1.580 * DL**3 * (DL + 4.) 
     3    + NF* (- 9.649 * DL1**2 - 406.171 * Y**2  - 32.218 * (1.-Y)
     4           - 5.976 * DL**2  - 1.60 * DL**3 )
       G2NSPA2 = 147.692 * DL1**2 + 2602.738 * Y**2 + 170.11 
     2           - 148.47 * DL    - 1.580 * DL**3 * (DL - 4.)
     3     + NF* ( 89.941 * DL1   - 218.482 * Y**2  - 9.623 
     4           - 0.910* DL**2   + 1.60 * DL**3 )
*
       IF (IMODN .EQ. 1) THEN
         G2NSPA = G2NSPA1
       ELSE IF (IMODN .EQ. 2) THEN
         G2NSPA = G2NSPA2
       ELSE
         G2NSPA = 0.5 * (G2NSPA1 + G2NSPA2)
       END IF
*
       G2NSPA = G2NSPA 
     1     + NF**2* (- 32.* Y*DL/(1.-Y) * (3.* DL + 10.) - 64.
     2               - (48.* DL**2 + 352.* DL + 384.) * (1.-Y) )/81.D0
       P2NSPA = -G2NSPA
*
       RETURN
       END
*
* ---------------------------------------------------------------------
*
* ..This is the regular piece of P^(2)-. 
*
       FUNCTION P2NSMA (Y, NF, IMODN)
*
       IMPLICIT REAL*8 (A-Z)
       INTEGER IMODN, NF
       DL  = LOG (Y)
       DL1 = LOG (1.-Y)
*
       G2NSMA1 =  157.387 * DL1**2 + 2741.42 * Y**2 + 490.43 * (1.-Y)
     2           - 67.0   * DL**2  - 10.005 * DL**3 - 1.432 * DL**4
     3    + NF *(- 17.989 * DL1**2 - 355.636 * Y**2 + 73.407 *(1.-Y)*DL1
     4           - 11.491 * DL**2  - 1.928 * DL**3 )
       G2NSMA2 = -115.099 * DL1**2 - 1581.05 * DL1  - 267.33  * (1.-Y)
     2           + 127.65 * DL**2  + 25.22 * DL**3  - 1.432 * DL**4
     3    + NF *(- 11.999 * DL1**2 - 397.546 * Y**2 - 41.949 * (1.-Y)
     4           + 1.477  * DL**2  + 0.538 * DL**3 )
*
       IF (IMODN .EQ. 1) THEN
         G2NSMA = G2NSMA1
       ELSE IF (IMODN .EQ. 2) THEN
         G2NSMA = G2NSMA2
       ELSE
         G2NSMA = 0.5 * (G2NSMA1 + G2NSMA2)
       END IF
*
       G2NSMA =  G2NSMA 
     3     + NF**2* (- 32.* Y*DL/(1.-Y) * (3.* DL + 10.) - 64.
     4               - (48.* DL**2 + 352.* DL + 384.) * (1.-Y) )/81.D0
       P2NSMA = -G2NSMA
*
       RETURN
       END
*
* ---------------------------------------------------------------------
*
* ..This is the singular NS piece for the `+' case. 
*
       FUNCTION P2NSPB (Y, NF, IMODN)
*
       IMPLICIT REAL*8 (A-Z)
       INTEGER IMODN, NF
       D1 = 1./(1.-Y)
*
       G2NSB1 = (-1183.762 + NF * 183.148) * D1
       G2NSB2 = (-1182.774 + NF * 183.931) * D1
*
       IF (IMODN .EQ. 1) THEN
         G2NSB = G2NSB1
       ELSE IF (IMODN .EQ. 2) THEN
         G2NSB = G2NSB2
       ELSE
         G2NSB = 0.5 * (G2NSB1 + G2NSB2)
       END IF
*
       G2NSB  =  G2NSB + NF**2 * 64./81.D0 * D1
       P2NSPB = -G2NSB
*
       RETURN
       END
*
* ---------------------------------------------------------------------
*
* ..This is the singular NS piece for the `-' case. 
*
       FUNCTION P2NSMB (Y, NF, IMODN)
*
       IMPLICIT REAL*8 (A-Z)
       DIMENSION CDM0(2), CDM1(2)
       INTEGER IMODN, NF
       DATA CDM0 / -1185.229, -1174.348 /
       DATA CDM1 /  184.765,  183.718 /
       D1 = 1./(1.-Y)
*
       G2NSB1 = (-1185.229 + NF * 184.765) * D1
       G2NSB2 = (-1174.348 + NF * 183.718) * D1
*
       IF (IMODN .EQ. 1) THEN
         G2NSB = G2NSB1
       ELSE IF (IMODN .EQ. 2) THEN
         G2NSB = G2NSB2
       ELSE
         G2NSB = 0.5 * (G2NSB1 + G2NSB2)
       END IF
*
       G2NSB  =  G2NSB + NF**2 * 64./81.D0 * D1
       P2NSMB = -G2NSB
*
       RETURN
       END
*
* ---------------------------------------------------------------------
*
* ..This is the 'local' NS+ piece piece for the `+' case.
*
       FUNCTION P2NSPC (Y, NF, IMODN)
*
       IMPLICIT REAL*8 (A-Z)
       INTEGER IMODN, NF
       DATA Z2, Z3 / 1.644934067 D0, 1.202056903 D0 /
       DL1 = LOG (1.-Y)
*
       G2NSC1 = (-1183.762 + NF* 183.148)* DL1 - 1347.032 + NF* 174.402
       G2NSC2 = (-1182.774 + NF* 183.931)* DL1 - 1351.088 + NF* 178.208
*
       IF (IMODN .EQ. 1) THEN
         G2NSC = G2NSC1
       ELSE IF (IMODN .EQ. 2) THEN
         G2NSC = G2NSC2
       ELSE
         G2NSC = 0.5 * (G2NSC1 + G2NSC2)
       END IF
*
       G2NSC  =  G2NSC
     1    + NF**2 * (16.* DL1 + (51. + 48.* Z3 - 80.* Z2)) * 4./81.D0
       P2NSPC = -G2NSC
*
       RETURN
       END
*
* ---------------------------------------------------------------------
*
* ..This is the 'local' NS+ piece piece for the `-' case.
*
       FUNCTION P2NSMC (Y, NF, IMODN)
*
       IMPLICIT REAL*8 (A-Z)
       INTEGER IMODN, NF
       DATA Z2, Z3 / 1.644934067 D0, 1.202056903 D0 /
       DL1 = LOG (1.-Y)
*
       G2NSC1 = (-1185.229 + NF* 184.765)* DL1 - 1365.458 + NF* 184.289
       G2NSC2 = (-1174.348 + NF* 183.718)* DL1 - 1286.799 + NF* 177.762
*
       IF (IMODN .EQ. 1) THEN
         G2NSC = G2NSC1
       ELSE IF (IMODN .EQ. 2) THEN
         G2NSC = G2NSC2
       ELSE
         G2NSC = 0.5 * (G2NSC1 + G2NSC2)
       END IF
*
       G2NSC  =  G2NSC
     1    + NF**2 * (16.* DL1 + (51. + 48.* Z3 - 80.* Z2)) * 4./81.D0
       P2NSMC = -G2NSC
*
       RETURN
       END
*
* ---------------------------------------------------------------------
*
* ..This is P^(2)S, the difference of P^(2)V and P^(2)-.
*
       FUNCTION P2NSSA (Y, NF, IMODN)
*
       IMPLICIT REAL*8 (A-Z)
       INTEGER IMODN, NF
       DL  = LOG (Y)
*
       G2NSSA1 = (  1441.57 * Y**2 - 12603.59 * Y + 15450.01) * (1.-Y)
     1            - 7876.93 * Y*DL**2 + 4260.29 * DL + 229.27 * DL**2
     2            - 4.4075 * DL**3
       G2NSSA2 = (  704.67 * Y**3 - 3310.32  * Y**2 - 2144.81  * Y
     1            + 244.68) * (1.-Y) - 4490.81 * Y**2 * DL
     2            - 42.875 * DL + 11.0165 * DL**3
*
       IF (IMODN .EQ. 1) THEN
         G2NSSA = G2NSSA1
       ELSE IF (IMODN .EQ. 2) THEN
         G2NSSA = G2NSSA2
       ELSE
         G2NSSA = 0.5 * (G2NSSA1 + G2NSSA2)
       END IF
*
       P2NSSA  = - NF * G2NSSA
*
       RETURN
       END
*
* =================================================================av==
