Files
HQQ/pphtt.f
Basil Bruhn 89acc8dd93 initial
Signed-off-by: Basil Bruhn <basil.bruhn@psi.ch>
2025-11-17 10:15:57 +01:00

2521 lines
102 KiB
Fortran

C**********************************************************************
C
C ***************
C * VERSION 1.2 *
C ***************
C
C--This program calculates the production cross section of Higgs
C bosons via qqbar,gg -> h,H,A + ttbar/bbbar at hadron colliders at LO
C QCD according to the results presented in
C
C Z. Kunszt, Nucl. Phys. B247 (1984) 339;
C J.F. Gunion, Phys. Lett. B261 (1991) 510;
C W.J. Marciano and F.E. Paige, Phys. Rev. Lett. 66 (1991) 2433;
C D.A. Dicus and S. Willenbrock, Phys. Rev. D39 (1989) 751;
C M. Spira, Fortschr. Phys. 46 (1998) 203.
C
C The program allows to calculate the total production cross sections
C for the neutral Higgs bosons of the SM and MSSM. The MSSM Higgs
C sector is implemented in the approximate two-loop RGE approach of
C
C M. Carena, J. Espinosa, M. Quiros and C.E.M. Wagner, Phys. Lett. B355
C (1995) 209.
C
C The program allows to use on-shell quark masses as well as MSbar
C quark masses at the scale of the invariant mass of the produced
C triplett of the Higgs boson and quarks for the Yukawa couplings. The
C parton densities are defined in the subroutine STRUC at the end of
C the program and can be changed there. The default is the use of the
C PDFLIB.
C
C--Description of the input file pphtt.in with default values:
C ===========================================================
C
C--Model: MSSM = 0 -> SM
C MSSM = 1 -> MSSM
C MSSM = 0
C
C--tan(beta) for MSSM
C TGBET = 3.D0
C
C--Higgs type: (SM: HIGGS = 2)
C HIGGS = 0 -> A
C HIGGS = 1 -> h
C HIGGS = 2 -> H
C HIGGS = 1
C
C--Higgs mass loop: MA = MA1 -> MA2 with NA equidistant points in total
C (SM: MH = MA)
C MA1 = 100.D0
C MA2 = 1000.D0
C NA = 10
C
C--Quark type (0=top, 1=bottom) and mass [GeV]:
C T=0/B=1 = 0
C MQ = 175.D0
C
C--Yukawa coupling: 0 = pole mass, 1 = running mass at scale Q = (M_H+2*M_Q)/2
C YUKAWA = 0
C
C--Energy of hadron collider in GeV
C ENERGY = 14000.D0
C
C--Choice pp (0) or ppbar (1) collider
C PP/PPBAR = 0
C
C--Cuts in the quark transverse momentum p_t and rapidity y:
C PTCUT = 0.D0
C YMIN = -100.D0
C YMAX = 100.D0
C
C--Scale loop: mu = xi * mu_0, xi = SCALE1 -> SCALE2 with NSCALE
C equidistant points in total
C mu_0 = Q = (M_H+2*M_Q)/2
C SCALE1 = 1.0D0
C SCALE2 = 1.0D0
C NSCALE = 1
C
C--Integration cut (should be left as it is)
C EPSILON = 1.D-8
C
C--Parton densities from PDFLIB (NPTYPE = 1). Here: CTEQ6L1
C NGROUP = -1
C NSET = 4
C Special cases: NGROUP = -1 --> CTEQ6
C -2 --> MRST2001
C
C--VEGAS: IPOINT = number of points
C ITERAT = number of iterations
C NPRN: 0 = no output of individual iterations
C 1 = full output of individual iterations
C 10 = table output of individual iterations
C IPOINT = 100000
C ITERAT = 5
C NPRN = 10
C
C--Heavy quark masses for running masses and alpha_s (top quark
C decoupled!!!)
C MC = 1.5D0
C MB = 5.D0
C MT = 175.D0
C
C--LAMBDA = LAMBDA_N0 for alpha_s in partonic cross sections
C N0 = 5
C LAMBDA = 0.181D0
C
C--Order of alpha_s: LOOP = 1 -> LO alpha_s
C LOOP = 2 -> NLO alpha_s
C LOOP = 1
C
C--Fermi constant, W mass, Z mass
C GF = 1.16639D-5
C MW = 80.41D0
C MZ = 91.187D0
C
C--MSSM parameters: in GeV
C MSQ = common squark mass
C MUR = right-handed up-type squark mass
C AT = trilinear coupling A_t
C AB = trilinear coupling A_b
C MU = Higgsino mass paramater mu
C MSQ = 1000.D0
C MUR = 1000.D0
C AT = 0.D0
C AB = 0.D0
C MU = 0.D0
C
C--LAMBDA = LAMBDA_N0 for alpha_s (NLO) used in MSSM Higgs sector
C N0 = 5
C LAMBDA = 0.202D0
C
C
C written by Michael Spira, e-mail: spira@desy.de.
C
C September 21, 1998
C=======================================================
C
PROGRAM PPHTT
C--PROGRAM FOR gg,qqbar -> h,H,A + QQbar AT HADRON COLLIDERS
PARAMETER(NIN=95,NOUT=96)
IMPLICIT DOUBLE PRECISION (A-H,O-Z)
COMMON/PARAM/S,AM,AMQ,SCALE,EPS,GF,GEVPB
COMMON/VMASSES/AMW,AMZ
COMMON/ALS/XLAMBDA,AMC,AMB,AMT,N0
COMMON/RESULT/RES,ERR,DCHI2,DUM
COMMON/PDFLIB/NGROUP,NSET
COMMON/FLAG/NHIGGS,LOOP
COMMON/COLLIDER/ICOLL
COMMON/CUTS/PTCUT,YMIN,YMAX
COMMON/YUKAWA/NYUKAWA,NQUARK
COMMON/MASSES/AMS0,AMC0,AMB0,AMT0
COMMON/STRANGE/AMSB
COMMON/HMASS/AMSM,AMA,AML,AMH,AMCH
COMMON/BREAK/AMSQ,AMUR,AU,AD,AMU,AM2,ASL
COMMON/COUP/GAT,GAB,GLT,GLB,GHT,GHB,GZAH,GZAL,
. GHHH,GLLL,GHLL,GLHH,GHAA,GLAA,GLVV,GHVV,
. GLPM,GHPM,B,A
EXTERNAL DSIGQQ,DSIGGG
PI = 4*DATAN(1.D0)
C--CONVERSION FACTOR 1/GEV**2 --> PB
GEVPB = 389379660.D0
OPEN(NIN,FILE='pphtt.in')
OPEN(NOUT,FILE='pphtt.out')
C--READING INPUT FILE
READ(NIN,101)IMODEL
READ(NIN,100)TGBET
READ(NIN,101)IHIGGS
IF(IMODEL.EQ.0) IHIGGS = 2
NHIGGS = 0
IF(IHIGGS.EQ.0) NHIGGS = 1
READ(NIN,100)AM1
READ(NIN,100)AM2
READ(NIN,101)NAM
READ(NIN,101)NQUARK
READ(NIN,100)AMQ
READ(NIN,101)NYUKAWA
READ(NIN,100)W
READ(NIN,101)NCOLL
C--NCOLL = 0 -> ICOLL = 1 -> PP COLLIDER
C--NCOLL = 1 -> ICOLL = -1 -> PPBAR COLLIDER
ICOLL=1
IF(NCOLL.EQ.1)ICOLL=-1
READ(NIN,100)PTCUT
READ(NIN,100)YMIN
READ(NIN,100)YMAX
READ(NIN,*)
READ(NIN,100)SCALE1
READ(NIN,100)SCALE2
READ(NIN,101)NSCALE
READ(NIN,100)EPS
READ(NIN,*)
READ(NIN,101)NGROUP
READ(NIN,101)NSET
READ(NIN,*)
READ(NIN,101)IPOINT
READ(NIN,101)ITERAT
READ(NIN,101)NPRN
READ(NIN,*)
READ(NIN,100)AMC0
READ(NIN,100)AMB0
READ(NIN,100)AMT0
READ(NIN,101)N00
READ(NIN,100)XLAMBDA0
READ(NIN,101)LOOP
READ(NIN,100)GF
READ(NIN,100)AMW
READ(NIN,100)AMZ
READ(NIN,*)
READ(NIN,100)AMSQ
READ(NIN,100)AMUR
READ(NIN,100)AU
READ(NIN,100)AD
READ(NIN,100)AMU
READ(NIN,101)N0S
READ(NIN,100)XLAMBDAS
C--Parameters not needed
ASL = 0
AMSB = 0.190D0
C--Parameters for COMMON block MASSES
AMS0 = AMSB
AMC = AMC0
AMB = AMB0
C--DECOUPLING TOP FROM RUNNING OF ALPHA_S
AMT = AMT0*1.D8
C--VEGAS PARAMETERS
ABSERR = 0.D0
IGRAPH = 0
C--INTIALIZING ALPHA_S
XLAMBDA = XLAMBDA0
N0 = N00
ACC = 1.D-10
CALL ALSINI(ACC)
C--INTIALIZING RANDOM NUMBER GENERATOR
CALL RSTART(12,34,56,78)
S = W**2
C--WRITING HEADER OF OUTPUT FILE
IF(NQUARK.EQ.0)THEN
IF(IHIGGS.EQ.0)THEN
WRITE(NOUT,*)' PP ---> ATT + X'
ELSEIF(IHIGGS.EQ.1)THEN
WRITE(NOUT,*)' PP ---> hTT + X'
ELSE
WRITE(NOUT,*)' PP ---> HTT + X'
ENDIF
ELSE
IF(IHIGGS.EQ.0)THEN
WRITE(NOUT,*)' PP ---> ABB + X'
ELSEIF(IHIGGS.EQ.1)THEN
WRITE(NOUT,*)' PP ---> hBB + X'
ELSE
WRITE(NOUT,*)' PP ---> HBB + X'
ENDIF
ENDIF
WRITE(NOUT,*)' ==============='
WRITE(NOUT,*)
WRITE(NOUT,*)' NGROUP = ',NGROUP,' NSET = ',NSET
WRITE(NOUT,*)' W = ',W,' GEV'
IF(NQUARK.EQ.0)THEN
WRITE(NOUT,*)' MT = ',AMQ,' GEV'
ELSE
WRITE(NOUT,*)' MB = ',AMQ,' GEV'
ENDIF
WRITE(NOUT,*)' LOOP = ',LOOP
IF(N0.EQ.4)THEN
WRITE(NOUT,*)' LAMBDA_4 = ',XLAMBDA
ELSE
WRITE(NOUT,*)' LAMBDA_5 = ',XLAMBDA
ENDIF
WRITE(NOUT,*)
. ' ----------------------------------------------------'
WRITE(NOUT,*)
DO 9999 I = 1,NAM
IF(NAM.EQ.1)THEN
AM = AM1
ELSE
AM = AM1 + (AM2-AM1)/(NAM-1.D0)*(I-1)
ENDIF
AMA = AM
COUPLING = 1.D0
C--Call MSSM couplings, if MSSM = 1
IF(IMODEL.EQ.1) THEN
C--INTIALIZING ALPHA_S FOR MSSM SUBROUTINE SUSYCP
XLAMBDA = XLAMBDAS
N0 = N0S
ACC = 1.D-10
CALL ALSINI(ACC)
CALL SUSYCP(TGBET)
C--INTIALIZING ALPHA_S FOR CROSS SECTION
XLAMBDA = XLAMBDA0
N0 = N00
ACC = 1.D-10
CALL ALSINI(ACC)
IF(NQUARK.EQ.0)THEN
IF(IHIGGS.EQ.0) THEN
AM = AMA
COUPLING = GAT**2
ELSEIF(IHIGGS.EQ.1) THEN
AM = AML
COUPLING = GLT**2
ELSE
AM = AMH
COUPLING = GHT**2
ENDIF
ELSE
IF(IHIGGS.EQ.0) THEN
AM = AMA
COUPLING = GAB**2
ELSEIF(IHIGGS.EQ.1) THEN
AM = AML
COUPLING = GLB**2
ELSE
AM = AMH
COUPLING = GHB**2
ENDIF
ENDIF
ENDIF
DO 9999 J = 1,NSCALE
IF(NSCALE.EQ.1)THEN
SCALE = SCALE1
ELSE
SCALE = SCALE1 + (SCALE2-SCALE1)/(NSCALE-1.D0)*(J-1)
ENDIF
C--INTEGRATION OF QQBAR INITIAL STATE
IDIM = 6
CALL VEGASN(DSIGQQ,ABSERR,IDIM,IPOINT,ITERAT,NPRN,IGRAPH)
SIGQ = RES*COUPLING
DSIGQ = ERR*COUPLING
C--INTEGRATION OF GG INITIAL STATE
CALL VEGASN(DSIGGG,ABSERR,IDIM,IPOINT,ITERAT,NPRN,IGRAPH)
SIGG = RES*COUPLING
DSIGG = ERR*COUPLING
SIGB = SIGQ + SIGG
DSIGB = DSQRT(DSIGQ**2 + DSIGG**2)
IF(IMODEL.EQ.0)THEN
WRITE(NOUT,*)' MH_SM = ',AM,' GEV'
ELSE
WRITE(NOUT,*)' MA_MSSM = ',AMA,' GEV'
WRITE(NOUT,*)' TG(BETA) = ',TGBET
IF(IHIGGS.EQ.0)THEN
WRITE(NOUT,*)' MA_MSSM = ',AMA,' GEV'
WRITE(NOUT,*)' G_Q^A**2 = ',COUPLING
ELSEIF(IHIGGS.EQ.1)THEN
WRITE(NOUT,*)' Mh_MSSM = ',AM,' GEV'
WRITE(NOUT,*)' G_Q^h**2 = ',COUPLING
ELSE
WRITE(NOUT,*)' MH_MSSM = ',AM,' GEV'
WRITE(NOUT,*)' G_Q^H**2 = ',COUPLING
ENDIF
ENDIF
WRITE(NOUT,*)' SCALE = ',SCALE
WRITE(NOUT,*)' SIGMA_QQ = (',SIGQ,' +- ',DSIGQ,') PB'
WRITE(NOUT,*)' SIGMA_GG = (',SIGG,' +- ',DSIGG,') PB'
WRITE(NOUT,*)' SIGMA = (',SIGB,' +- ',DSIGB,') PB'
WRITE(NOUT,*)
9999 CONTINUE
100 FORMAT(10X,G30.20)
101 FORMAT(10X,I30)
STOP
END
DOUBLE PRECISION FUNCTION DSIGQQ(XX)
C--FUNCTION FOR QQBAR INITIAL STATE
IMPLICIT DOUBLE PRECISION (A-H,O-Z)
COMMON/PARAM/S,AM,AMQ,SCALE,EPS,GF,GEVPB
DIMENSION XX(6),Y(6)
DO I=1,6
Y(I) = EPS + (1-2*EPS)*XX(I)
ENDDO
DSIGQQ = GEVPB*QMAT(Y)*(1-2*EPS)**6
RETURN
END
DOUBLE PRECISION FUNCTION DSIGGG(XX)
C--FUNCTION FOR GG INITIAL STATE
IMPLICIT DOUBLE PRECISION (A-H,O-Z)
COMMON/PARAM/S,AM,AMQ,SCALE,EPS,GF,GEVPB
DIMENSION XX(6),Y(6)
DO I=1,6
Y(I) = EPS + (1-2*EPS)*XX(I)
ENDDO
DSIGGG = GEVPB*GMAT(Y)*(1-2*EPS)**6
RETURN
END
DOUBLE PRECISION FUNCTION QMAT(Y)
C--QQBAR CROSS SECTION AND KINEMATICS
IMPLICIT DOUBLE PRECISION (A-H,O-Z)
DIMENSION PDF(-6:6), Y(6)
COMMON/PARAM/S,AM,AMQ,SCALE,EPS,GF,GEVPB
COMMON/FLAG/NHIGGS,LOOP
COMMON/CUTS/PTCUT,YMIN,YMAX
COMMON/YUKAWA/NYUKAWA,NQUARK
QMAT = 0
PI = 4*DATAN(1.D0)
TAUH = (AM+2*AMQ)**2/S
TAU = TAUH + (1-TAUH)*Y(1)
Y1 = -DLOG(TAU)*Y(2)
X1 = DEXP(-Y1)
X2 = TAU/X1
SHAT = X1*X2*S
c>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
c QQ = SCALE*DSQRT(SHAT)
QQ = SCALE*(AM/2+AMQ)
c>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
C--RESCALING YUKAWA COUPLING FOR RUNNING QUARK MASS
RAT = 1
IF(NYUKAWA.NE.0)THEN
IF(NQUARK.EQ.0)THEN
RAT = (RUNM(QQ,6)/AMQ)**2
ELSE
RAT = (RUNM(QQ,5)/AMQ)**2
ENDIF
ENDIF
ALP = ALPHAS(QQ,LOOP)/PI
COUP = ALP**2/1152.D0*DSQRT(2.D0)*GF*AMQ**2 * RAT
XMU1 = AMQ**2/SHAT
XMU2 = AMQ**2/SHAT
XMU3 = AM**2/SHAT
XX3P = 1+XMU3-(DSQRT(XMU1)+DSQRT(XMU2))**2
XX3M = 2*DSQRT(XMU3)
XX3 = XX3M + (XX3P-XX3M)*Y(3)
AA = 1-XX3+XMU3
BB = (2-XX3)*(AA+XMU1-XMU2)
DELTA = (XX3**2-4*XMU3)*(AA-(DSQRT(XMU1)+DSQRT(XMU2))**2)
. *(AA-(DSQRT(XMU1)-DSQRT(XMU2))**2)
XX1M = (BB-DSQRT(DELTA))/2/AA
XX1P = (BB+DSQRT(DELTA))/2/AA
XX1 = XX1M + (XX1P-XX1M)*Y(4)
XX2 = 2-XX1-XX3
FAC0 = (XX3P-XX3M)*(XX1P-XX1M)
CTH = -1 + 2*Y(5)
CHI = 2*PI*Y(6)
FAC = -(1-TAUH)*DLOG(TAU)*2*2*PI * FAC0
BET1 = DSQRT(1-4*AMQ**2/XX1**2/SHAT)
BET2 = DSQRT(1-4*AMQ**2/XX2**2/SHAT)
CT12 = 2/XX1/XX2/BET1/BET2*(1-XX1-XX2+XX1*XX2/2+2*XMU1-XMU3)
ST12 = DSQRT(1-CT12**2)
CCHI = DCOS(CHI)
STH = DSQRT(1-CTH**2)
X12 = SHAT
X13 = -SHAT/2*XX1*(1-BET1*CTH)
X23 = -SHAT/2*XX1*(1+BET1*CTH)
X14 = -SHAT/2*XX2*(1-BET2*(ST12*CCHI*STH+CT12*CTH))
X24 = -SHAT/2*XX2*(1+BET2*(ST12*CCHI*STH+CT12*CTH))
IF(NHIGGS.NE.1)THEN
QMAT = COUP*DHQQ(X12,X13,X14,X23,X24)*DLUMQQ(TAU,X1,QQ)
. * FAC
ELSE
QMAT = COUP*DAQQ(X12,X13,X14,X23,X24)*DLUMQQ(TAU,X1,QQ)
. * FAC
ENDIF
PTCUT2 = PTCUT**2
PT1SQ = BET1**2*XX1**2*SHAT/4*STH**2
PT2SQ = BET2**2*XX2**2*SHAT/4*(1-(ST12*CCHI*STH+CT12*CTH)**2)
RAP1 = DLOG((1+BET1*CTH)/(1-BET1*CTH))/2
RAP2 = DLOG((1+BET2*(ST12*CCHI*STH+CT12*CTH))/
. (1-BET2*(ST12*CCHI*STH+CT12*CTH)))/2
IF(PT1SQ.LE.PTCUT2.OR.PT2SQ.LE.PTCUT2) QMAT = 0
IF(RAP1.LE.YMIN.OR.RAP1.GE.YMAX) QMAT = 0
IF(RAP2.LE.YMIN.OR.RAP2.GE.YMAX) QMAT = 0
RETURN
END
DOUBLE PRECISION FUNCTION GMAT(Y)
C--GG CROSS SECTION AND KINEMATICS
IMPLICIT DOUBLE PRECISION (A-H,O-Z)
DIMENSION PDF(-6:6), Y(6)
COMMON/PARAM/S,AM,AMQ,SCALE,EPS,GF,GEVPB
COMMON/FLAG/NHIGGS,LOOP
COMMON/CUTS/PTCUT,YMIN,YMAX
COMMON/YUKAWA/NYUKAWA,NQUARK
GMAT = 0
PI = 4*DATAN(1.D0)
TAUH = (AM+2*AMQ)**2/S
TAU = TAUH + (1-TAUH)*Y(1)
Y1 = -DLOG(TAU)*Y(2)
X1 = DEXP(-Y1)
X2 = TAU/X1
SHAT = X1*X2*S
c>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
c QQ = SCALE*DSQRT(SHAT)
QQ = SCALE*(AM/2+AMQ)
c>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
C--RESCALING YUKAWA COUPLING FOR RUNNING QUARK MASS
RAT = 1
IF(NYUKAWA.NE.0)THEN
IF(NQUARK.EQ.0)THEN
RAT = (RUNM(QQ,6)/AMQ)**2
ELSE
RAT = (RUNM(QQ,5)/AMQ)**2
ENDIF
ENDIF
ALP = ALPHAS(QQ,LOOP)/PI
COUP = ALP**2/6144.D0*DSQRT(2.D0)*GF*AMQ**2 * RAT
XMU1 = AMQ**2/SHAT
XMU2 = AMQ**2/SHAT
XMU3 = AM**2/SHAT
XX3P = 1+XMU3-(DSQRT(XMU1)+DSQRT(XMU2))**2
XX3M = 2*DSQRT(XMU3)
XX3 = XX3M + (XX3P-XX3M)*Y(3)
AA = 1-XX3+XMU3
BB = (2-XX3)*(AA+XMU1-XMU2)
DELTA = (XX3**2-4*XMU3)*(AA-(DSQRT(XMU1)+DSQRT(XMU2))**2)
. *(AA-(DSQRT(XMU1)-DSQRT(XMU2))**2)
XX1M = (BB-DSQRT(DELTA))/2/AA
XX1P = (BB+DSQRT(DELTA))/2/AA
XX1 = XX1M + (XX1P-XX1M)*Y(4)
XX2 = 2-XX1-XX3
FAC0 = (XX3P-XX3M)*(XX1P-XX1M)
CTH = -1 + 2*Y(5)
CHI = 2*PI*Y(6)
FAC = -(1-TAUH)*DLOG(TAU)*2*2*PI * FAC0
BET1 = DSQRT(1-4*AMQ**2/XX1**2/SHAT)
BET2 = DSQRT(1-4*AMQ**2/XX2**2/SHAT)
CT12 = 2/XX1/XX2/BET1/BET2*(1-XX1-XX2+XX1*XX2/2+2*XMU1-XMU3)
ST12 = DSQRT(1-CT12**2)
CCHI = DCOS(CHI)
STH = DSQRT(1-CTH**2)
X12 = SHAT
X13 = -SHAT/2*XX1*(1-BET1*CTH)
X23 = -SHAT/2*XX1*(1+BET1*CTH)
X14 = -SHAT/2*XX2*(1-BET2*(ST12*CCHI*STH+CT12*CTH))
X24 = -SHAT/2*XX2*(1+BET2*(ST12*CCHI*STH+CT12*CTH))
IF(NHIGGS.NE.1)THEN
GMAT = COUP*DHGG(X12,X13,X14,X23,X24)*DLUMGG(TAU,X1,QQ)
. * FAC
ELSE
GMAT = COUP*DAGG(X12,X13,X14,X23,X24)*DLUMGG(TAU,X1,QQ)
. * FAC
ENDIF
PTCUT2 = PTCUT**2
PT1SQ = BET1**2*XX1**2*SHAT/4*STH**2
PT2SQ = BET2**2*XX2**2*SHAT/4*(1-(ST12*CCHI*STH+CT12*CTH)**2)
RAP1 = DLOG((1+BET1*CTH)/(1-BET1*CTH))/2
RAP2 = DLOG((1+BET2*(ST12*CCHI*STH+CT12*CTH))/
. (1-BET2*(ST12*CCHI*STH+CT12*CTH)))/2
IF(PT1SQ.LE.PTCUT2.OR.PT2SQ.LE.PTCUT2) GMAT = 0
IF(RAP1.LE.YMIN.OR.RAP1.GE.YMAX) GMAT = 0
IF(RAP2.LE.YMIN.OR.RAP2.GE.YMAX) GMAT = 0
RETURN
END
DOUBLE PRECISION FUNCTION DHQQ(X12,X13,X14,X23,X24)
C--QQBAR MATRIX ELEMENT FOR SCALAR HIGGS BOSONS
IMPLICIT DOUBLE PRECISION (A-Z)
COMMON/PARAM/S,AM,AMQ,SCALE,EPS,GF,GEVPB
N12 = 1/X12
N13 = 1/X13
N14 = 1/X14
N23 = 1/X23
N24 = 1/X24
N35 = 1/(X12+X14+X24)
N45 = 1/(X12+X13+X23)
M = AMQ
MH = AM
RES =
+ - 96*M**2*MH**2*N35*N45*X12 - 16*M**2*MH**2*N35**2*X12 - 16*
+ M**2*MH**2*N45**2*X12 - 64*M**2*N35*N45*X12*X13 - 64*M**2*N35*
+ N45*X12*X14 - 128*M**2*N35*N45*X13*X14 + 64*M**2*N35*X12 + 64*
+ M**2*N35**2*X12*X14 + 32*M**2*N35**2*X12**2 + 64*M**2*N35**2*
+ X14**2 + 64*M**2*N45*X12 + 64*M**2*N45**2*X12*X13 + 32*M**2*
+ N45**2*X12**2 + 64*M**2*N45**2*X13**2 + 128*M**4*N35*N45*X12 +
+ 64*M**4*N35**2*X12 + 64*M**4*N45**2*X12 + 16*MH**2*N35*N45*X12*
+ X13 + 16*MH**2*N35*N45*X12*X14 + 16*MH**2*N35*N45*X12**2 + 32*
+ MH**2*N35*N45*X13*X14 - 16*MH**2*N35*X12 - 16*MH**2*N35**2*X12*
+ X14 - 8*MH**2*N35**2*X12**2 - 16*MH**2*N35**2*X14**2 - 16*MH**2*
+ N45*X12 - 16*MH**2*N45**2*X12*X13 - 8*MH**2*N45**2*X12**2 - 16*
+ MH**2*N45**2*X13**2 + 16*MH**4*N35*N45*X12 + 32*N35*N45*X12*X13*
+ X14 + 16*N35*N45*X12*X13**2 + 16*N35*N45*X12*X14**2 + 32*N35*N45
+ *X12**2*X13 + 32*N35*N45*X12**2*X14 + 16*N35*N45*X12**3 - 8*N35*
+ X12*X13
RES = RES - 16*N35*X12*X14 + 8*N35*X12*X23 - 8*N35*X12**2 - 16*
+ N45*X12*X13 - 8*N45*X12*X14 + 8*N45*X12*X24 - 8*N45*X12**2 + 16*
+ X12
DHQQ = RES*N12**2
RETURN
END
DOUBLE PRECISION FUNCTION DAQQ(X12,X13,X14,X23,X24)
C--QQBAR MATRIX ELEMENT FOR PSEUDOSCALAR HIGGS BOSONS
IMPLICIT DOUBLE PRECISION (A-Z)
COMMON/PARAM/S,AM,AMQ,SCALE,EPS,GF,GEVPB
N12 = 1/X12
N13 = 1/X13
N14 = 1/X14
N23 = 1/X23
N24 = 1/X24
N35 = 1/(X12+X14+X24)
N45 = 1/(X12+X13+X23)
M = AMQ
MH = AM
RES =
+ - 32*M**2*MH**2*N35*N45*X12 - 16*M**2*MH**2*N35**2*X12 - 16*
+ M**2*MH**2*N45**2*X12 + 16*MH**2*N35*N45*X12*X13 + 16*MH**2*N35*
+ N45*X12*X14 + 16*MH**2*N35*N45*X12**2 + 32*MH**2*N35*N45*X13*X14
+ - 16*MH**2*N35*X12 - 16*MH**2*N35**2*X12*X14 - 8*MH**2*N35**2*
+ X12**2 - 16*MH**2*N35**2*X14**2 - 16*MH**2*N45*X12 - 16*MH**2*
+ N45**2*X12*X13 - 8*MH**2*N45**2*X12**2 - 16*MH**2*N45**2*X13**2
+ + 16*MH**4*N35*N45*X12 + 32*N35*N45*X12*X13*X14 + 16*N35*N45*
+ X12*X13**2 + 16*N35*N45*X12*X14**2 + 32*N35*N45*X12**2*X13 + 32*
+ N35*N45*X12**2*X14 + 16*N35*N45*X12**3 - 8*N35*X12*X13 - 16*N35*
+ X12*X14 + 8*N35*X12*X23 - 8*N35*X12**2 - 16*N45*X12*X13 - 8*N45*
+ X12*X14 + 8*N45*X12*X24 - 8*N45*X12**2 + 16*X12
DAQQ = RES*N12**2
RETURN
END
DOUBLE PRECISION FUNCTION DHGG(X12,X13,X14,X23,X24)
C--GG MATRIX ELEMENT FOR SCALAR HIGGS BOSONS
IMPLICIT DOUBLE PRECISION (A-Z)
COMMON/PARAM/S,AM,AMQ,SCALE,EPS,GF,GEVPB
N12 = 1/X12
N13 = 1/X13
N14 = 1/X14
N23 = 1/X23
N24 = 1/X24
N35 = 1/(X12+X14+X24)
N45 = 1/(X12+X13+X23)
M = AMQ
MH = AM
RES =
+ 108*M**2*MH**2*N12*N13*N24*N35*X14 - 8*M**2*MH**2*N12*N13*N24*
+ N35*X23 - 72*M**2*MH**2*N12*N13*N24*N45*X14 - 324*M**2*MH**2*N12
+ *N13*N24 + 64*M**2*MH**2*N12*N13*N24**2*N35*X14*X23 - 64*M**2*
+ MH**2*N12*N13*N24**2*X23 + 72*M**2*MH**2*N12*N13*N35*N45*X14 -
+ 108*M**2*MH**2*N12*N13*N35 - 72*M**2*MH**2*N12*N13*N45 - 72*M**2
+ *MH**2*N12*N14*N23*N35*X13 + 108*M**2*MH**2*N12*N14*N23*N45*X13
+ - 8*M**2*MH**2*N12*N14*N23*N45*X24 - 324*M**2*MH**2*N12*N14*N23
+ + 64*M**2*MH**2*N12*N14*N23**2*N45*X13*X24 - 64*M**2*MH**2*N12*
+ N14*N23**2*X24 + 72*M**2*MH**2*N12*N14*N35*N45*X13 - 72*M**2*
+ MH**2*N12*N14*N35 - 108*M**2*MH**2*N12*N14*N45 - 432*M**2*MH**2*
+ N12*N23*N35*N45*X13 - 72*M**2*MH**2*N12*N23*N35*N45*X14 + 216*
+ M**2*MH**2*N12*N23*N35 + 92*M**2*MH**2*N12*N23*N45 - 92*M**2*
+ MH**2*N12*N23*N45**2*X13 - 72*M**2*MH**2*N12*N24*N35*N45*X13 -
+ 432*M**2*MH**2*N12*N24*N35*N45*X14 + 92*M**2*MH**2*N12*N24*N35
+ - 92*M**2*MH**2*N12*N24*N35**2*X14
RES = RES + 216*M**2*MH**2*N12*N24*N45 + 864*M**2*MH**2*N12*N35*
+ N45 + 52*M**2*MH**2*N12*N35**2 + 52*M**2*MH**2*N12*N45**2 - 40*
+ M**2*MH**2*N13*N14*N23*N24*X12 - 4*M**2*MH**2*N13*N14*N23*N45*
+ X12 - 8*M**2*MH**2*N13*N14*N23*N45*X24 - 52*M**2*MH**2*N13*N14*
+ N23 - 4*M**2*MH**2*N13*N14*N24*N35*X12 - 8*M**2*MH**2*N13*N14*
+ N24*N35*X23 - 52*M**2*MH**2*N13*N14*N24 + 40*M**2*MH**2*N13*N14*
+ N35*N45*X12 - 52*M**2*MH**2*N13*N14*N35 - 52*M**2*MH**2*N13*N14*
+ N45 - 4*M**2*MH**2*N13*N23*N24*N45*X12 - 8*M**2*MH**2*N13*N23*
+ N24*N45*X14 - 52*M**2*MH**2*N13*N23*N24 - 24*M**2*MH**2*N13*N23*
+ N45 + 8*M**2*MH**2*N13*N23*N45**2*X12 + 96*M**2*MH**2*N13*N24*
+ N35*N45*X12 + 32*M**2*MH**2*N13*N24*N35*N45*X14 + 392*M**2*MH**2
+ *N13*N24*N35 + 156*M**2*MH**2*N13*N24*N45 + 96*M**2*MH**2*N13*
+ N24**2*N35*X12 + 96*M**2*MH**2*N13*N24**2*N35*X14 + 64*M**2*
+ MH**2*N13*N24**2*N35*X23 - 32*M**2*MH**2*N13*N24**2 + 472*M**2*
+ MH**2*N13*N35*N45 + 72*M**2*MH**2*N13*N45**2 + 64*M**2*MH**2*
+ N13**2*N24
RES = RES + 64*M**2*MH**2*N13**2*N45 - 4*M**2*MH**2*N14*N23*N24*
+ N35*X12 - 8*M**2*MH**2*N14*N23*N24*N35*X13 - 52*M**2*MH**2*N14*
+ N23*N24 + 96*M**2*MH**2*N14*N23*N35*N45*X12 + 32*M**2*MH**2*N14*
+ N23*N35*N45*X13 + 156*M**2*MH**2*N14*N23*N35 + 392*M**2*MH**2*
+ N14*N23*N45 + 96*M**2*MH**2*N14*N23**2*N45*X12 + 96*M**2*MH**2*
+ N14*N23**2*N45*X13 + 64*M**2*MH**2*N14*N23**2*N45*X24 - 32*M**2*
+ MH**2*N14*N23**2 - 24*M**2*MH**2*N14*N24*N35 + 8*M**2*MH**2*N14*
+ N24*N35**2*X12 + 472*M**2*MH**2*N14*N35*N45 + 72*M**2*MH**2*N14*
+ N35**2 + 64*M**2*MH**2*N14**2*N23 + 64*M**2*MH**2*N14**2*N35 -
+ 64*M**2*MH**2*N23*N24*N35*N45*X12 - 52*M**2*MH**2*N23*N24*N35*
+ N45*X13 - 52*M**2*MH**2*N23*N24*N35*N45*X14 - 116*M**2*MH**2*N23
+ *N35*N45 - 52*M**2*MH**2*N23*N45**2 + 96*M**2*MH**2*N23**2*N45
+ - 32*M**2*MH**2*N23**2*N45**2*X12 - 32*M**2*MH**2*N23**2*N45**2
+ *X13 - 116*M**2*MH**2*N24*N35*N45 - 52*M**2*MH**2*N24*N35**2 +
+ 96*M**2*MH**2*N24**2*N35 - 32*M**2*MH**2*N24**2*N35**2*X12 - 32*
+ M**2*MH**2*N24**2*N35**2*X14
RES = RES + 32*M**2*MH**4*N13*N14*N23*N24 + 8*M**2*MH**4*N13*N14*
+ N23*N45 + 8*M**2*MH**4*N13*N14*N24*N35 + 32*M**2*MH**4*N13*N14*
+ N35*N45 + 8*M**2*MH**4*N13*N23*N24*N45 - 256*M**2*MH**4*N13*N24*
+ N35*N45 - 64*M**2*MH**4*N13*N24**2*N35 - 64*M**2*MH**4*N13**2*
+ N24*N45 + 8*M**2*MH**4*N14*N23*N24*N35 - 256*M**2*MH**4*N14*N23*
+ N35*N45 - 64*M**2*MH**4*N14*N23**2*N45 - 64*M**2*MH**4*N14**2*
+ N23*N35 + 32*M**2*MH**4*N23*N24*N35*N45 - 34*M**2*N12*N13*N24*
+ N35*X14*X23 + 6*M**2*N12*N13*N24*N35*X14**2 + 138*M**2*N12*N13*
+ N24*X14 + 202*M**2*N12*N13*N24*X23 - 32*M**2*N12*N13*N24**2*N35*
+ X14**2*X23 + 32*M**2*N12*N13*N24**2*X14*X23 - 138*M**2*N12*N13*
+ N35*X14 + 86*M**2*N12*N13*N35*X23 + 288*M**2*N12*N13 - 34*M**2*
+ N12*N14*N23*N45*X13*X24 + 6*M**2*N12*N14*N23*N45*X13**2 + 138*
+ M**2*N12*N14*N23*X13 + 202*M**2*N12*N14*N23*X24 - 32*M**2*N12*
+ N14*N23**2*N45*X13**2*X24 + 32*M**2*N12*N14*N23**2*X13*X24 - 138
+ *M**2*N12*N14*N45*X13 + 86*M**2*N12*N14*N45*X24 + 288*M**2*N12*
+ N14
RES = RES + 72*M**2*N12*N23*N35*N45*X13*X14 + 52*M**2*N12*N23*N35
+ *N45*X13**2 + 92*M**2*N12*N23*N35*X13 - 72*M**2*N12*N23*N35*X14
+ + 52*M**2*N12*N23*N45*X13 - 14*M**2*N12*N23*N45*X14 - 74*M**2*
+ N12*N23*N45*X24 + 14*M**2*N12*N23*N45**2*X13*X14 + 42*M**2*N12*
+ N23*N45**2*X13*X24 + 46*M**2*N12*N23*N45**2*X13**2 + 190*M**2*
+ N12*N23 - 96*M**2*N12*N23**2*N45*X13*X24 + 32*M**2*N12*N23**2*
+ N45**2*X13**2*X24 + 64*M**2*N12*N23**2*X24 + 72*M**2*N12*N24*N35
+ *N45*X13*X14 + 52*M**2*N12*N24*N35*N45*X14**2 - 14*M**2*N12*N24*
+ N35*X13 + 52*M**2*N12*N24*N35*X14 - 74*M**2*N12*N24*N35*X23 + 14
+ *M**2*N12*N24*N35**2*X13*X14 + 42*M**2*N12*N24*N35**2*X14*X23 +
+ 46*M**2*N12*N24*N35**2*X14**2 - 72*M**2*N12*N24*N45*X13 + 92*
+ M**2*N12*N24*N45*X14 + 190*M**2*N12*N24 - 96*M**2*N12*N24**2*N35
+ *X14*X23 + 32*M**2*N12*N24**2*N35**2*X14**2*X23 + 64*M**2*N12*
+ N24**2*X23 + 556*M**2*N12*N35*N45*X13 + 556*M**2*N12*N35*N45*X14
+ - 622*M**2*N12*N35 + 14*M**2*N12*N35**2*X13 - 242*M**2*N12*
+ N35**2*X14
RES = RES + 10*M**2*N12*N35**2*X23 - 622*M**2*N12*N45 - 242*M**2*
+ N12*N45**2*X13 + 14*M**2*N12*N45**2*X14 + 10*M**2*N12*N45**2*X24
+ + 56*M**2*N12**2*N13*N24*N35*X14**2*X23 - 56*M**2*N12**2*N13*
+ N24*X14*X23 + 56*M**2*N12**2*N13*N35*X14*X23 + 56*M**2*N12**2*
+ N14*N23*N45*X13**2*X24 - 56*M**2*N12**2*N14*N23*X13*X24 + 56*
+ M**2*N12**2*N14*N45*X13*X24 - 144*M**2*N12**2*N23*N35*N45*X13**2
+ *X14 + 144*M**2*N12**2*N23*N35*X13*X14 + 144*M**2*N12**2*N23*N45
+ *X13**2 - 144*M**2*N12**2*N23*X13 - 144*M**2*N12**2*N24*N35*N45*
+ X13*X14**2 + 144*M**2*N12**2*N24*N35*X14**2 + 144*M**2*N12**2*
+ N24*N45*X13*X14 - 144*M**2*N12**2*N24*X14 + 288*M**2*N12**2*N35*
+ N45*X13*X14 + 144*M**2*N12**2*N35*X14 - 288*M**2*N12**2*N35**2*
+ X14**2 + 144*M**2*N12**2*N45*X13 - 288*M**2*N12**2*N45**2*X13**2
+ + 12*M**2*N13*N14*N23*N24*X12**2 + 32*M**2*N13*N14*N23*X12 + 20
+ *M**2*N13*N14*N23*X24 + 32*M**2*N13*N14*N24*X12 + 20*M**2*N13*
+ N14*N24*X23 + 12*M**2*N13*N14*N35*N45*X12**2 - 12*M**2*N13*N14*
+ N35*X12
RES = RES + 20*M**2*N13*N14*N35*X23 - 12*M**2*N13*N14*N45*X12 +
+ 20*M**2*N13*N14*N45*X24 + 80*M**2*N13*N14 + 32*M**2*N13*N23*N24*
+ X12 + 20*M**2*N13*N23*N24*X14 + 8*M**2*N13*N23*N45*X12 + 8*M**2*
+ N13*N23*N45*X14 + 8*M**2*N13*N23*N45*X24 - 4*M**2*N13*N23*N45**2
+ *X12*X14 - 4*M**2*N13*N23*N45**2*X12*X24 - 4*M**2*N13*N23*N45**2
+ *X12**2 + 32*M**2*N13*N23 - 264*M**2*N13*N24*N35*X12 - 162*M**2*
+ N13*N24*N35*X14 - 122*M**2*N13*N24*N35*X23 - 96*M**2*N13*N24*N45
+ *X12 - 32*M**2*N13*N24*N45*X14 + 220*M**2*N13*N24 - 64*M**2*N13*
+ N24**2*N35*X12*X14 - 32*M**2*N13*N24**2*N35*X12*X23 - 32*M**2*
+ N13*N24**2*N35*X12**2 - 64*M**2*N13*N24**2*N35*X14*X23 - 32*M**2
+ *N13*N24**2*N35*X14**2 + 32*M**2*N13*N24**2*X23 + 128*M**2*N13*
+ N35*N45*X12 + 180*M**2*N13*N35*N45*X14 - 368*M**2*N13*N35 - 236*
+ M**2*N13*N45 - 68*M**2*N13*N45**2*X12 - 4*M**2*N13*N45**2*X14 -
+ 4*M**2*N13*N45**2*X24 - 32*M**2*N13**2*N24*X12 - 32*M**2*N13**2*
+ N24*X23 - 32*M**2*N13**2*N45*X24 - 64*M**2*N13**2 + 32*M**2*N14*
+ N23*N24*X12
RES = RES + 20*M**2*N14*N23*N24*X13 - 96*M**2*N14*N23*N35*X12 -
+ 32*M**2*N14*N23*N35*X13 - 264*M**2*N14*N23*N45*X12 - 162*M**2*
+ N14*N23*N45*X13 - 122*M**2*N14*N23*N45*X24 + 220*M**2*N14*N23 -
+ 64*M**2*N14*N23**2*N45*X12*X13 - 32*M**2*N14*N23**2*N45*X12*X24
+ - 32*M**2*N14*N23**2*N45*X12**2 - 64*M**2*N14*N23**2*N45*X13*
+ X24 - 32*M**2*N14*N23**2*N45*X13**2 + 32*M**2*N14*N23**2*X24 + 8
+ *M**2*N14*N24*N35*X12 + 8*M**2*N14*N24*N35*X13 + 8*M**2*N14*N24*
+ N35*X23 - 4*M**2*N14*N24*N35**2*X12*X13 - 4*M**2*N14*N24*N35**2*
+ X12*X23 - 4*M**2*N14*N24*N35**2*X12**2 + 32*M**2*N14*N24 + 128*
+ M**2*N14*N35*N45*X12 + 180*M**2*N14*N35*N45*X13 - 236*M**2*N14*
+ N35 - 68*M**2*N14*N35**2*X12 - 4*M**2*N14*N35**2*X13 - 4*M**2*
+ N14*N35**2*X23 - 368*M**2*N14*N45 - 32*M**2*N14**2*N23*X12 - 32*
+ M**2*N14**2*N23*X24 - 32*M**2*N14**2*N35*X23 - 64*M**2*N14**2 +
+ 56*M**2*N23*N24*N35*N45*X12*X13 + 56*M**2*N23*N24*N35*N45*X12*
+ X14 + 36*M**2*N23*N24*N35*N45*X12**2 + 48*M**2*N23*N24*N35*N45*
+ X13*X14
RES = RES + 20*M**2*N23*N24*N35*N45*X13**2 + 20*M**2*N23*N24*N35*
+ N45*X14**2 - 8*M**2*N23*N24*N35*X12 - 8*M**2*N23*N24*N35*X14 - 8
+ *M**2*N23*N24*N45*X12 - 8*M**2*N23*N24*N45*X13 + 48*M**2*N23*N24
+ + 68*M**2*N23*N35*N45*X12 + 184*M**2*N23*N35*N45*X13 + 56*M**2*
+ N23*N35*N45*X14 - 180*M**2*N23*N35 - 290*M**2*N23*N45 + 56*M**2*
+ N23*N45**2*X12 + 170*M**2*N23*N45**2*X13 + 10*M**2*N23*N45**2*
+ X14 + 38*M**2*N23*N45**2*X24 - 128*M**2*N23**2*N45*X12 - 128*
+ M**2*N23**2*N45*X13 - 32*M**2*N23**2*N45*X14 - 96*M**2*N23**2*
+ N45*X24 + 64*M**2*N23**2*N45**2*X12*X13 + 32*M**2*N23**2*N45**2*
+ X12*X24 + 32*M**2*N23**2*N45**2*X12**2 + 64*M**2*N23**2*N45**2*
+ X13*X24 + 32*M**2*N23**2*N45**2*X13**2 + 32*M**2*N23**2 + 68*
+ M**2*N24*N35*N45*X12 + 56*M**2*N24*N35*N45*X13 + 184*M**2*N24*
+ N35*N45*X14 - 290*M**2*N24*N35 + 56*M**2*N24*N35**2*X12 + 10*
+ M**2*N24*N35**2*X13 + 170*M**2*N24*N35**2*X14 + 38*M**2*N24*
+ N35**2*X23 - 180*M**2*N24*N45 - 128*M**2*N24**2*N35*X12 - 32*
+ M**2*N24**2*N35*X13
RES = RES - 128*M**2*N24**2*N35*X14 - 96*M**2*N24**2*N35*X23 + 64
+ *M**2*N24**2*N35**2*X12*X14 + 32*M**2*N24**2*N35**2*X12*X23 + 32
+ *M**2*N24**2*N35**2*X12**2 + 64*M**2*N24**2*N35**2*X14*X23 + 32*
+ M**2*N24**2*N35**2*X14**2 + 32*M**2*N24**2 + 680*M**2*N35*N45 -
+ 180*M**2*N35**2 - 180*M**2*N45**2 - 80*M**4*MH**2*N13*N14*N23*
+ N24 - 48*M**4*MH**2*N13*N14*N23*N45 - 48*M**4*MH**2*N13*N14*N24*
+ N35 - 80*M**4*MH**2*N13*N14*N35*N45 - 48*M**4*MH**2*N13*N23*N24*
+ N45 - 16*M**4*MH**2*N13*N23*N45**2 + 640*M**4*MH**2*N13*N24*N35*
+ N45 + 384*M**4*MH**2*N13*N24**2*N35 + 384*M**4*MH**2*N13**2*N24*
+ N45 + 64*M**4*MH**2*N13**2*N24**2 + 64*M**4*MH**2*N13**2*N45**2
+ - 48*M**4*MH**2*N14*N23*N24*N35 + 640*M**4*MH**2*N14*N23*N35*
+ N45 + 384*M**4*MH**2*N14*N23**2*N45 - 16*M**4*MH**2*N14*N24*
+ N35**2 + 384*M**4*MH**2*N14**2*N23*N35 + 64*M**4*MH**2*N14**2*
+ N23**2 + 64*M**4*MH**2*N14**2*N35**2 - 80*M**4*MH**2*N23*N24*N35
+ *N45 + 64*M**4*MH**2*N23**2*N45**2 + 64*M**4*MH**2*N24**2*N35**2
+ - 144*M**4*N12*N13*N24*N35*X14
RES = RES + 32*M**4*N12*N13*N24*N35*X23 + 288*M**4*N12*N13*N24*
+ N45*X14 + 432*M**4*N12*N13*N24 - 256*M**4*N12*N13*N24**2*N35*X14
+ *X23 + 256*M**4*N12*N13*N24**2*X23 - 288*M**4*N12*N13*N35*N45*
+ X14 + 144*M**4*N12*N13*N35 + 288*M**4*N12*N13*N45 + 288*M**4*N12
+ *N14*N23*N35*X13 - 144*M**4*N12*N14*N23*N45*X13 + 32*M**4*N12*
+ N14*N23*N45*X24 + 432*M**4*N12*N14*N23 - 256*M**4*N12*N14*N23**2
+ *N45*X13*X24 + 256*M**4*N12*N14*N23**2*X24 - 288*M**4*N12*N14*
+ N35*N45*X13 + 288*M**4*N12*N14*N35 + 144*M**4*N12*N14*N45 + 576*
+ M**4*N12*N23*N35*N45*X13 + 288*M**4*N12*N23*N35*N45*X14 - 288*
+ M**4*N12*N23*N35 - 368*M**4*N12*N23*N45 + 368*M**4*N12*N23*
+ N45**2*X13 + 288*M**4*N12*N24*N35*N45*X13 + 576*M**4*N12*N24*N35
+ *N45*X14 - 368*M**4*N12*N24*N35 + 368*M**4*N12*N24*N35**2*X14 -
+ 288*M**4*N12*N24*N45 - 1152*M**4*N12*N35*N45 - 208*M**4*N12*
+ N35**2 - 208*M**4*N12*N45**2 + 64*M**4*N13*N14*N23*N24*X12 + 32*
+ M**4*N13*N14*N23*N45*X12 + 32*M**4*N13*N14*N23*N45*X24 + 64*M**4
+ *N13*N14*N23
RES = RES + 32*M**4*N13*N14*N24*N35*X12 + 32*M**4*N13*N14*N24*N35
+ *X23 + 64*M**4*N13*N14*N24 - 64*M**4*N13*N14*N35*N45*X12 + 96*
+ M**4*N13*N14*N35 + 96*M**4*N13*N14*N45 + 32*M**4*N13*N23*N24*N45
+ *X12 + 32*M**4*N13*N23*N24*N45*X14 + 64*M**4*N13*N23*N24 + 64*
+ M**4*N13*N23*N45 - 624*M**4*N13*N24*N35 - 224*M**4*N13*N24*N45
+ - 256*M**4*N13*N24**2*N35*X12 - 256*M**4*N13*N24**2*N35*X14 -
+ 256*M**4*N13*N24**2*N35*X23 - 608*M**4*N13*N35*N45 - 256*M**4*
+ N13*N45**2 - 256*M**4*N13**2*N24 - 256*M**4*N13**2*N45 + 32*M**4
+ *N14*N23*N24*N35*X12 + 32*M**4*N14*N23*N24*N35*X13 + 64*M**4*N14
+ *N23*N24 - 224*M**4*N14*N23*N35 - 624*M**4*N14*N23*N45 - 256*
+ M**4*N14*N23**2*N45*X12 - 256*M**4*N14*N23**2*N45*X13 - 256*M**4
+ *N14*N23**2*N45*X24 + 64*M**4*N14*N24*N35 - 608*M**4*N14*N35*N45
+ - 256*M**4*N14*N35**2 - 256*M**4*N14**2*N23 - 256*M**4*N14**2*
+ N35 + 64*M**4*N23*N24*N35*N45*X12 + 64*M**4*N23*N24*N35*N45*X13
+ + 64*M**4*N23*N24*N35*N45*X14 + 32*M**4*N23*N24*N35 + 32*M**4*
+ N23*N24*N45
RES = RES + 320*M**4*N23*N35*N45 + 112*M**4*N23*N45**2 - 256*M**4
+ *N23**2*N45 + 320*M**4*N24*N35*N45 + 112*M**4*N24*N35**2 - 256*
+ M**4*N24**2*N35 + 64*M**6*N13*N14*N23*N24 + 64*M**6*N13*N14*N23*
+ N45 + 64*M**6*N13*N14*N24*N35 + 64*M**6*N13*N14*N35*N45 + 64*
+ M**6*N13*N23*N24*N45 + 64*M**6*N13*N23*N45**2 - 512*M**6*N13*N24
+ *N35*N45 - 512*M**6*N13*N24**2*N35 - 512*M**6*N13**2*N24*N45 -
+ 256*M**6*N13**2*N24**2 - 256*M**6*N13**2*N45**2 + 64*M**6*N14*
+ N23*N24*N35 - 512*M**6*N14*N23*N35*N45 - 512*M**6*N14*N23**2*N45
+ + 64*M**6*N14*N24*N35**2 - 512*M**6*N14**2*N23*N35 - 256*M**6*
+ N14**2*N23**2 - 256*M**6*N14**2*N35**2 + 64*M**6*N23*N24*N35*N45
+ - 256*M**6*N23**2*N45**2 - 256*M**6*N24**2*N35**2 - MH**2*N12*
+ N13*N24*N35*X14*X23 - 5*MH**2*N12*N13*N24*N35*X14**2 - 31*MH**2*
+ N12*N13*N24*X14 - 49*MH**2*N12*N13*N24*X23 + 31*MH**2*N12*N13*
+ N35*X14 - 23*MH**2*N12*N13*N35*X23 - 72*MH**2*N12*N13 - MH**2*
+ N12*N14*N23*N45*X13*X24 - 5*MH**2*N12*N14*N23*N45*X13**2 - 31*
+ MH**2*N12*N14*N23*X13
RES = RES - 49*MH**2*N12*N14*N23*X24 + 31*MH**2*N12*N14*N45*X13
+ - 23*MH**2*N12*N14*N45*X24 - 72*MH**2*N12*N14 - 18*MH**2*N12*
+ N23*N35*N45*X13*X14 - 40*MH**2*N12*N23*N35*N45*X13**2 + 4*MH**2*
+ N12*N23*N35*X13 + 18*MH**2*N12*N23*N35*X14 - 41*MH**2*N12*N23*
+ N45*X13 + 23*MH**2*N12*N23*N45**2*X13**2 - 54*MH**2*N12*N23 - 18
+ *MH**2*N12*N24*N35*N45*X13*X14 - 40*MH**2*N12*N24*N35*N45*X14**2
+ - 41*MH**2*N12*N24*N35*X14 + 23*MH**2*N12*N24*N35**2*X14**2 +
+ 18*MH**2*N12*N24*N45*X13 + 4*MH**2*N12*N24*N45*X14 - 54*MH**2*
+ N12*N24 - 166*MH**2*N12*N35*N45*X13 - 166*MH**2*N12*N35*N45*X14
+ + 162*MH**2*N12*N35 + 95*MH**2*N12*N35**2*X14 + 162*MH**2*N12*
+ N45 + 95*MH**2*N12*N45**2*X13 - 14*MH**2*N12**2*N13*N24*N35*
+ X14**2*X23 + 14*MH**2*N12**2*N13*N24*X14*X23 - 14*MH**2*N12**2*
+ N13*N35*X14*X23 - 14*MH**2*N12**2*N14*N23*N45*X13**2*X24 + 14*
+ MH**2*N12**2*N14*N23*X13*X24 - 14*MH**2*N12**2*N14*N45*X13*X24
+ + 36*MH**2*N12**2*N23*N35*N45*X13**2*X14 - 36*MH**2*N12**2*N23*
+ N35*X13*X14
RES = RES - 36*MH**2*N12**2*N23*N45*X13**2 + 36*MH**2*N12**2*N23*
+ X13 + 36*MH**2*N12**2*N24*N35*N45*X13*X14**2 - 36*MH**2*N12**2*
+ N24*N35*X14**2 - 36*MH**2*N12**2*N24*N45*X13*X14 + 36*MH**2*
+ N12**2*N24*X14 - 72*MH**2*N12**2*N35*N45*X13*X14 - 36*MH**2*
+ N12**2*N35*X14 + 72*MH**2*N12**2*N35**2*X14**2 - 36*MH**2*N12**2
+ *N45*X13 + 72*MH**2*N12**2*N45**2*X13**2 - 6*MH**2*N13*N14*N23*
+ N24*X12**2 + 2*MH**2*N13*N14*N23*N45*X12*X24 + 2*MH**2*N13*N14*
+ N23*N45*X12**2 - 14*MH**2*N13*N14*N23*X12 - 8*MH**2*N13*N14*N23*
+ X24 + 2*MH**2*N13*N14*N24*N35*X12*X23 + 2*MH**2*N13*N14*N24*N35*
+ X12**2 - 14*MH**2*N13*N14*N24*X12 - 8*MH**2*N13*N14*N24*X23 - 6*
+ MH**2*N13*N14*N35*N45*X12**2 + 8*MH**2*N13*N14*N35*X12 - 4*MH**2
+ *N13*N14*N35*X23 + 8*MH**2*N13*N14*N45*X12 - 4*MH**2*N13*N14*N45
+ *X24 - 24*MH**2*N13*N14 + 2*MH**2*N13*N23*N24*N45*X12*X14 + 2*
+ MH**2*N13*N23*N24*N45*X12**2 - 14*MH**2*N13*N23*N24*X12 - 8*
+ MH**2*N13*N23*N24*X14 + 8*MH**2*N13*N23*N45*X12 - 2*MH**2*N13*
+ N23*N45**2*X12**2
RES = RES - 18*MH**2*N13*N23 + 47*MH**2*N13*N24*N35*X12 + 24*
+ MH**2*N13*N24*N35*X14 + 31*MH**2*N13*N24*N35*X23 + 18*MH**2*N13*
+ N24*N45*X12 + 18*MH**2*N13*N24*N45*X14 - 77*MH**2*N13*N24 - 60*
+ MH**2*N13*N35*N45*X12 - 54*MH**2*N13*N35*N45*X14 + 105*MH**2*N13
+ *N35 + 60*MH**2*N13*N45 + 14*MH**2*N13*N45**2*X12 + 2*MH**2*N14*
+ N23*N24*N35*X12*X13 + 2*MH**2*N14*N23*N24*N35*X12**2 - 14*MH**2*
+ N14*N23*N24*X12 - 8*MH**2*N14*N23*N24*X13 + 18*MH**2*N14*N23*N35
+ *X12 + 18*MH**2*N14*N23*N35*X13 + 47*MH**2*N14*N23*N45*X12 + 24*
+ MH**2*N14*N23*N45*X13 + 31*MH**2*N14*N23*N45*X24 - 77*MH**2*N14*
+ N23 + 8*MH**2*N14*N24*N35*X12 - 2*MH**2*N14*N24*N35**2*X12**2 -
+ 18*MH**2*N14*N24 - 60*MH**2*N14*N35*N45*X12 - 54*MH**2*N14*N35*
+ N45*X13 + 60*MH**2*N14*N35 + 14*MH**2*N14*N35**2*X12 + 105*MH**2
+ *N14*N45 - 22*MH**2*N23*N24*N35*N45*X12*X13 - 22*MH**2*N23*N24*
+ N35*N45*X12*X14 - 14*MH**2*N23*N24*N35*N45*X12**2 - 20*MH**2*N23
+ *N24*N35*N45*X13*X14 - 8*MH**2*N23*N24*N35*N45*X13**2 - 8*MH**2*
+ N23*N24*N35*N45*X14**2
RES = RES + 10*MH**2*N23*N24*N35*X12 + 4*MH**2*N23*N24*N35*X13 +
+ 8*MH**2*N23*N24*N35*X14 + 10*MH**2*N23*N24*N45*X12 + 8*MH**2*N23
+ *N24*N45*X13 + 4*MH**2*N23*N24*N45*X14 - 20*MH**2*N23*N24 - 14*
+ MH**2*N23*N35*N45*X12 - 54*MH**2*N23*N35*N45*X13 - 8*MH**2*N23*
+ N35*N45*X14 + 42*MH**2*N23*N35 + 35*MH**2*N23*N45 + 3*MH**2*N23*
+ N45**2*X12 + 12*MH**2*N23*N45**2*X13 - 14*MH**2*N24*N35*N45*X12
+ - 8*MH**2*N24*N35*N45*X13 - 54*MH**2*N24*N35*N45*X14 + 35*MH**2
+ *N24*N35 + 3*MH**2*N24*N35**2*X12 + 12*MH**2*N24*N35**2*X14 + 42
+ *MH**2*N24*N45 - 232*MH**2*N35*N45 + 57*MH**2*N35**2 + 57*MH**2*
+ N45**2 - 18*MH**4*N12*N13*N24*N35*X14 + 54*MH**4*N12*N13*N24 +
+ 18*MH**4*N12*N13*N35 - 18*MH**4*N12*N14*N23*N45*X13 + 54*MH**4*
+ N12*N14*N23 + 18*MH**4*N12*N14*N45 + 72*MH**4*N12*N23*N35*N45*
+ X13 - 36*MH**4*N12*N23*N35 + 72*MH**4*N12*N24*N35*N45*X14 - 36*
+ MH**4*N12*N24*N45 - 144*MH**4*N12*N35*N45 + 8*MH**4*N13*N14*N23*
+ N24*X12 - 2*MH**4*N13*N14*N23*N45*X12 + 10*MH**4*N13*N14*N23 - 2
+ *MH**4*N13*N14*N24*N35*X12
RES = RES + 10*MH**4*N13*N14*N24 - 8*MH**4*N13*N14*N35*N45*X12 +
+ 6*MH**4*N13*N14*N35 + 6*MH**4*N13*N14*N45 - 2*MH**4*N13*N23*N24*
+ N45*X12 + 10*MH**4*N13*N23*N24 - 16*MH**4*N13*N24*N35*N45*X12 -
+ 16*MH**4*N13*N24*N35*N45*X14 - 52*MH**4*N13*N24*N35 - 18*MH**4*
+ N13*N24*N45 - 88*MH**4*N13*N35*N45 - 2*MH**4*N14*N23*N24*N35*X12
+ + 10*MH**4*N14*N23*N24 - 16*MH**4*N14*N23*N35*N45*X12 - 16*
+ MH**4*N14*N23*N35*N45*X13 - 18*MH**4*N14*N23*N35 - 52*MH**4*N14*
+ N23*N45 - 88*MH**4*N14*N35*N45 + 12*MH**4*N23*N24*N35*N45*X12 +
+ 10*MH**4*N23*N24*N35*N45*X13 + 10*MH**4*N23*N24*N35*N45*X14 - 4*
+ MH**4*N23*N24*N35 - 4*MH**4*N23*N24*N45 + 10*MH**4*N23*N35*N45
+ + 10*MH**4*N24*N35*N45 - 4*MH**6*N13*N14*N23*N24 - 4*MH**6*N13*
+ N14*N35*N45 + 32*MH**6*N13*N24*N35*N45 + 32*MH**6*N14*N23*N35*
+ N45 - 4*MH**6*N23*N24*N35*N45 - 2*N12*N13*N24*N35*X14*X23**2 +
+ 28*N12*N13*N24*N35*X14**2*X23 + 7*N12*N13*N24*N35*X14**3 + 15*
+ N12*N13*N24*X14*X23 + 11*N12*N13*N24*X14**2 + 27*N12*N13*N24*
+ X23**2
RES = RES - 15*N12*N13*N35*X14*X23 + 25*N12*N13*N35*X14**2 + 9*
+ N12*N13*N35*X23**2 + 18*N12*N13*X14 + 72*N12*N13*X23 + 36*N12*
+ N13*X24 - 2*N12*N14*N23*N45*X13*X24**2 + 28*N12*N14*N23*N45*
+ X13**2*X24 + 7*N12*N14*N23*N45*X13**3 + 15*N12*N14*N23*X13*X24
+ + 11*N12*N14*N23*X13**2 + 27*N12*N14*N23*X24**2 - 15*N12*N14*
+ N45*X13*X24 + 25*N12*N14*N45*X13**2 + 9*N12*N14*N45*X24**2 + 18*
+ N12*N14*X13 + 36*N12*N14*X23 + 72*N12*N14*X24 - 16*N12*N23*N35*
+ N45*X13**2*X14 + 2*N12*N23*N35*N45*X13**3 + 34*N12*N23*N35*X13*
+ X14 + 16*N12*N23*N35*X13**2 + 9*N12*N23*N45*X13*X14 + 45*N12*N23
+ *N45*X13*X24 + 98*N12*N23*N45*X13**2 + 9*N12*N23*N45**2*X13**2*
+ X14 - 35*N12*N23*N45**2*X13**2*X24 - 23*N12*N23*N45**2*X13**3 -
+ 21*N12*N23*X13 + 18*N12*N23*X14 + 33*N12*N23*X24 - 16*N12*N24*
+ N35*N45*X13*X14**2 + 2*N12*N24*N35*N45*X14**3 + 9*N12*N24*N35*
+ X13*X14 + 45*N12*N24*N35*X14*X23 + 98*N12*N24*N35*X14**2 + 9*N12
+ *N24*N35**2*X13*X14**2 - 35*N12*N24*N35**2*X14**2*X23 - 23*N12*
+ N24*N35**2*X14**3
RES = RES + 34*N12*N24*N45*X13*X14 + 16*N12*N24*N45*X14**2 + 18*
+ N12*N24*X13 - 21*N12*N24*X14 + 33*N12*N24*X23 - 140*N12*N35*N45*
+ X13*X14 - 52*N12*N35*N45*X13**2 - 52*N12*N35*N45*X14**2 + 36*N12
+ *N35*X13 + 129*N12*N35*X14 - 33*N12*N35*X23 + 9*N12*N35**2*X13*
+ X14 - 12*N12*N35**2*X14*X23 - 23*N12*N35**2*X14**2 + 129*N12*N45
+ *X13 + 36*N12*N45*X14 - 33*N12*N45*X24 + 9*N12*N45**2*X13*X14 -
+ 12*N12*N45**2*X13*X24 - 23*N12*N45**2*X13**2 + 7*N12**2*N13*N24*
+ N35*X14**2*X23**2 + 7*N12**2*N13*N24*N35*X14**3*X23 - 7*N12**2*
+ N13*N24*X14*X23**2 - 7*N12**2*N13*N24*X14**2*X23 + 7*N12**2*N13*
+ N35*X14*X23**2 + 7*N12**2*N13*N35*X14**2*X23 + 7*N12**2*N14*N23*
+ N45*X13**2*X24**2 + 7*N12**2*N14*N23*N45*X13**3*X24 - 7*N12**2*
+ N14*N23*X13*X24**2 - 7*N12**2*N14*N23*X13**2*X24 + 7*N12**2*N14*
+ N45*X13*X24**2 + 7*N12**2*N14*N45*X13**2*X24 - 18*N12**2*N23*N35
+ *N45*X13**2*X14**2 - 18*N12**2*N23*N35*N45*X13**3*X14 + 18*
+ N12**2*N23*N35*X13*X14**2 + 18*N12**2*N23*N35*X13**2*X14 + 18*
+ N12**2*N23*N45*X13**2*X14
RES = RES + 53*N12**2*N23*N45*X13**2*X24 + 18*N12**2*N23*N45*
+ X13**3 - 23*N12**2*N23*N45**2*X13**3*X24 - 18*N12**2*N23*X13*X14
+ - 30*N12**2*N23*X13*X24 - 18*N12**2*N23*X13**2 - 18*N12**2*N24*
+ N35*N45*X13*X14**3 - 18*N12**2*N24*N35*N45*X13**2*X14**2 + 18*
+ N12**2*N24*N35*X13*X14**2 + 53*N12**2*N24*N35*X14**2*X23 + 18*
+ N12**2*N24*N35*X14**3 - 23*N12**2*N24*N35**2*X14**3*X23 + 18*
+ N12**2*N24*N45*X13*X14**2 + 18*N12**2*N24*N45*X13**2*X14 - 18*
+ N12**2*N24*X13*X14 - 30*N12**2*N24*X14*X23 - 18*N12**2*N24*
+ X14**2 - 36*N12**2*N35*N45*X13*X14**2 - 36*N12**2*N35*N45*X13**2
+ *X14 + 18*N12**2*N35*X13*X14 + 30*N12**2*N35*X14*X23 + 18*N12**2
+ *N35*X14**2 - 23*N12**2*N35**2*X14**2*X23 + 18*N12**2*N45*X13*
+ X14 + 30*N12**2*N45*X13*X24 + 18*N12**2*N45*X13**2 - 23*N12**2*
+ N45**2*X13**2*X24 + 2*N13*N14*N23*N24*X12**3 + 6*N13*N14*N23*X12
+ *X24 + 6*N13*N14*N23*X12**2 + 2*N13*N14*N23*X24**2 + 6*N13*N14*
+ N24*X12*X23 + 6*N13*N14*N24*X12**2 + 2*N13*N14*N24*X23**2 - 2*
+ N13*N14*N35*N45*X12**3
RES = RES - 2*N13*N14*N35*X12*X23 + 2*N13*N14*N35*X12**2 + 2*N13*
+ N14*N35*X23**2 - 2*N13*N14*N45*X12*X24 + 2*N13*N14*N45*X12**2 +
+ 2*N13*N14*N45*X24**2 + 12*N13*N14*X12 + 12*N13*N14*X23 + 12*N13*
+ N14*X24 + 6*N13*N23*N24*X12*X14 + 6*N13*N23*N24*X12**2 + 2*N13*
+ N23*N24*X14**2 - 4*N13*N23*N45*X12*X14 - 4*N13*N23*N45*X12*X24
+ - 4*N13*N23*N45*X12**2 + 2*N13*N23*N45**2*X12**2*X14 + 2*N13*
+ N23*N45**2*X12**2*X24 + 2*N13*N23*N45**2*X12**3 + 14*N13*N23*X12
+ + 8*N13*N23*X14 + 8*N13*N23*X24 + 5*N13*N24*N35*X12*X14 - 18*
+ N13*N24*N35*X12*X23 - 9*N13*N24*N35*X12**2 + 3*N13*N24*N35*X14*
+ X23 + 21*N13*N24*N35*X14**2 - 9*N13*N24*N35*X23**2 + 39*N13*N24*
+ X12 + 28*N13*N24*X14 + 60*N13*N24*X23 - 38*N13*N35*N45*X12*X14
+ - 22*N13*N35*N45*X12**2 - 18*N13*N35*N45*X14**2 + 13*N13*N35*
+ X12 + 52*N13*N35*X14 - 40*N13*N35*X23 + 18*N13*N45*X12 + 16*N13*
+ N45*X14 - 20*N13*N45*X24 + 2*N13*N45**2*X12*X14 + 2*N13*N45**2*
+ X12*X24 + 2*N13*N45**2*X12**2 + 28*N13 + 6*N14*N23*N24*X12*X13
+ + 6*N14*N23*N24*X12**2
RES = RES + 2*N14*N23*N24*X13**2 + 5*N14*N23*N45*X12*X13 - 18*N14
+ *N23*N45*X12*X24 - 9*N14*N23*N45*X12**2 + 3*N14*N23*N45*X13*X24
+ + 21*N14*N23*N45*X13**2 - 9*N14*N23*N45*X24**2 + 39*N14*N23*X12
+ + 28*N14*N23*X13 + 60*N14*N23*X24 - 4*N14*N24*N35*X12*X13 - 4*
+ N14*N24*N35*X12*X23 - 4*N14*N24*N35*X12**2 + 2*N14*N24*N35**2*
+ X12**2*X13 + 2*N14*N24*N35**2*X12**2*X23 + 2*N14*N24*N35**2*
+ X12**3 + 14*N14*N24*X12 + 8*N14*N24*X13 + 8*N14*N24*X23 - 38*N14
+ *N35*N45*X12*X13 - 22*N14*N35*N45*X12**2 - 18*N14*N35*N45*X13**2
+ + 18*N14*N35*X12 + 16*N14*N35*X13 - 20*N14*N35*X23 + 2*N14*
+ N35**2*X12*X13 + 2*N14*N35**2*X12*X23 + 2*N14*N35**2*X12**2 + 13
+ *N14*N45*X12 + 52*N14*N45*X13 - 40*N14*N45*X24 + 28*N14 + 24*N23
+ *N24*N35*N45*X12*X13*X14 + 10*N23*N24*N35*N45*X12*X13**2 + 10*
+ N23*N24*N35*N45*X12*X14**2 + 14*N23*N24*N35*N45*X12**2*X13 + 14*
+ N23*N24*N35*N45*X12**2*X14 + 6*N23*N24*N35*N45*X12**3 + 10*N23*
+ N24*N35*N45*X13*X14**2 + 10*N23*N24*N35*N45*X13**2*X14 + 2*N23*
+ N24*N35*N45*X13**3
RES = RES + 2*N23*N24*N35*N45*X14**3 - 4*N23*N24*N35*X12*X13 - 8*
+ N23*N24*N35*X12*X14 - 4*N23*N24*N35*X12**2 - 4*N23*N24*N35*X13*
+ X14 - 4*N23*N24*N35*X14**2 - 8*N23*N24*N45*X12*X13 - 4*N23*N24*
+ N45*X12*X14 - 4*N23*N24*N45*X12**2 - 4*N23*N24*N45*X13*X14 - 4*
+ N23*N24*N45*X13**2 + 16*N23*N24*X12 + 10*N23*N24*X13 + 10*N23*
+ N24*X14 + 14*N23*N35*N45*X12*X13 + 8*N23*N35*N45*X12*X14 + 6*N23
+ *N35*N45*X12**2 + 10*N23*N35*N45*X13*X14 + 10*N23*N35*N45*X13**2
+ + 2*N23*N35*N45*X14**2 - 4*N23*N35*X12 + 12*N23*N35*X13 - 4*N23
+ *N35*X14 + 6*N23*N45*X12 + 74*N23*N45*X13 - 13*N23*N45*X14 + 4*
+ N23*N45*X24 - 15*N23*N45**2*X12*X13 + 13*N23*N45**2*X12*X14 -
+ N23*N45**2*X12*X24 - N23*N45**2*X12**2 + 20*N23*N45**2*X13*X14
+ - 15*N23*N45**2*X13*X24 - 35*N23*N45**2*X13**2 + 43*N23 + 8*N24
+ *N35*N45*X12*X13 + 14*N24*N35*N45*X12*X14 + 6*N24*N35*N45*X12**2
+ + 10*N24*N35*N45*X13*X14 + 2*N24*N35*N45*X13**2 + 10*N24*N35*
+ N45*X14**2 + 6*N24*N35*X12 - 13*N24*N35*X13 + 74*N24*N35*X14 + 4
+ *N24*N35*X23
RES = RES + 13*N24*N35**2*X12*X13 - 15*N24*N35**2*X12*X14 - N24*
+ N35**2*X12*X23 - N24*N35**2*X12**2 + 20*N24*N35**2*X13*X14 - 15*
+ N24*N35**2*X14*X23 - 35*N24*N35**2*X14**2 - 4*N24*N45*X12 - 4*
+ N24*N45*X13 + 12*N24*N45*X14 + 43*N24 - 108*N35*N45*X12 - 138*
+ N35*N45*X13 - 138*N35*N45*X14 + 95*N35 - 3*N35**2*X12 + 11*
+ N35**2*X13 - 12*N35**2*X14 - 3*N35**2*X23 + 95*N45 - 3*N45**2*
+ X12 - 12*N45**2*X13 + 11*N45**2*X14 - 3*N45**2*X24
DHGG = RES
RETURN
END
DOUBLE PRECISION FUNCTION DAGG(X12,X13,X14,X23,X24)
C--GG MATRIX ELEMENT FOR PSEUDOSCALAR HIGGS BOSONS
IMPLICIT DOUBLE PRECISION (A-Z)
COMMON/PARAM/S,AM,AMQ,SCALE,EPS,GF,GEVPB
N12 = 1/X12
N13 = 1/X13
N14 = 1/X14
N23 = 1/X23
N24 = 1/X24
N35 = 1/(X12+X14+X24)
N45 = 1/(X12+X13+X23)
M = AMQ
MH = AM
RES =
+ 36*M**2*MH**2*N12*N13*N24*N35*X14 - 8*M**2*MH**2*N12*N13*N24*N35
+ *X23 - 72*M**2*MH**2*N12*N13*N24*N45*X14 - 108*M**2*MH**2*N12*
+ N13*N24 + 64*M**2*MH**2*N12*N13*N24**2*N35*X14*X23 - 64*M**2*
+ MH**2*N12*N13*N24**2*X23 + 72*M**2*MH**2*N12*N13*N35*N45*X14 -
+ 36*M**2*MH**2*N12*N13*N35 - 72*M**2*MH**2*N12*N13*N45 - 72*M**2*
+ MH**2*N12*N14*N23*N35*X13 + 36*M**2*MH**2*N12*N14*N23*N45*X13 -
+ 8*M**2*MH**2*N12*N14*N23*N45*X24 - 108*M**2*MH**2*N12*N14*N23 +
+ 64*M**2*MH**2*N12*N14*N23**2*N45*X13*X24 - 64*M**2*MH**2*N12*N14
+ *N23**2*X24 + 72*M**2*MH**2*N12*N14*N35*N45*X13 - 72*M**2*MH**2*
+ N12*N14*N35 - 36*M**2*MH**2*N12*N14*N45 - 144*M**2*MH**2*N12*N23
+ *N35*N45*X13 - 72*M**2*MH**2*N12*N23*N35*N45*X14 + 72*M**2*MH**2
+ *N12*N23*N35 + 92*M**2*MH**2*N12*N23*N45 - 92*M**2*MH**2*N12*N23
+ *N45**2*X13 - 72*M**2*MH**2*N12*N24*N35*N45*X13 - 144*M**2*MH**2
+ *N12*N24*N35*N45*X14 + 92*M**2*MH**2*N12*N24*N35 - 92*M**2*MH**2
+ *N12*N24*N35**2*X14
RES = RES + 72*M**2*MH**2*N12*N24*N45 + 288*M**2*MH**2*N12*N35*
+ N45 + 52*M**2*MH**2*N12*N35**2 + 52*M**2*MH**2*N12*N45**2 - 16*
+ M**2*MH**2*N13*N14*N23*N24*X12 - 4*M**2*MH**2*N13*N14*N23*N45*
+ X12 - 8*M**2*MH**2*N13*N14*N23*N45*X24 - 20*M**2*MH**2*N13*N14*
+ N23 - 4*M**2*MH**2*N13*N14*N24*N35*X12 - 8*M**2*MH**2*N13*N14*
+ N24*N35*X23 - 20*M**2*MH**2*N13*N14*N24 + 16*M**2*MH**2*N13*N14*
+ N35*N45*X12 - 20*M**2*MH**2*N13*N14*N35 - 20*M**2*MH**2*N13*N14*
+ N45 - 4*M**2*MH**2*N13*N23*N24*N45*X12 - 8*M**2*MH**2*N13*N23*
+ N24*N45*X14 - 20*M**2*MH**2*N13*N23*N24 - 24*M**2*MH**2*N13*N23*
+ N45 + 8*M**2*MH**2*N13*N23*N45**2*X12 + 32*M**2*MH**2*N13*N24*
+ N35*N45*X12 + 32*M**2*MH**2*N13*N24*N35*N45*X14 + 192*M**2*MH**2
+ *N13*N24*N35 + 28*M**2*MH**2*N13*N24*N45 + 96*M**2*MH**2*N13*
+ N24**2*N35*X12 + 96*M**2*MH**2*N13*N24**2*N35*X14 + 64*M**2*
+ MH**2*N13*N24**2*N35*X23 - 32*M**2*MH**2*N13*N24**2 + 184*M**2*
+ MH**2*N13*N35*N45 + 72*M**2*MH**2*N13*N45**2 + 64*M**2*MH**2*
+ N13**2*N24
RES = RES + 64*M**2*MH**2*N13**2*N45 - 4*M**2*MH**2*N14*N23*N24*
+ N35*X12 - 8*M**2*MH**2*N14*N23*N24*N35*X13 - 20*M**2*MH**2*N14*
+ N23*N24 + 32*M**2*MH**2*N14*N23*N35*N45*X12 + 32*M**2*MH**2*N14*
+ N23*N35*N45*X13 + 28*M**2*MH**2*N14*N23*N35 + 192*M**2*MH**2*N14
+ *N23*N45 + 96*M**2*MH**2*N14*N23**2*N45*X12 + 96*M**2*MH**2*N14*
+ N23**2*N45*X13 + 64*M**2*MH**2*N14*N23**2*N45*X24 - 32*M**2*
+ MH**2*N14*N23**2 - 24*M**2*MH**2*N14*N24*N35 + 8*M**2*MH**2*N14*
+ N24*N35**2*X12 + 184*M**2*MH**2*N14*N35*N45 + 72*M**2*MH**2*N14*
+ N35**2 + 64*M**2*MH**2*N14**2*N23 + 64*M**2*MH**2*N14**2*N35 -
+ 24*M**2*MH**2*N23*N24*N35*N45*X12 - 20*M**2*MH**2*N23*N24*N35*
+ N45*X13 - 20*M**2*MH**2*N23*N24*N35*N45*X14 - 84*M**2*MH**2*N23*
+ N35*N45 - 52*M**2*MH**2*N23*N45**2 + 96*M**2*MH**2*N23**2*N45 -
+ 32*M**2*MH**2*N23**2*N45**2*X12 - 32*M**2*MH**2*N23**2*N45**2*
+ X13 - 84*M**2*MH**2*N24*N35*N45 - 52*M**2*MH**2*N24*N35**2 + 96*
+ M**2*MH**2*N24**2*N35 - 32*M**2*MH**2*N24**2*N35**2*X12 - 32*
+ M**2*MH**2*N24**2*N35**2*X14
RES = RES + 16*M**2*MH**4*N13*N14*N23*N24 + 8*M**2*MH**4*N13*N14*
+ N23*N45 + 8*M**2*MH**4*N13*N14*N24*N35 + 16*M**2*MH**4*N13*N14*
+ N35*N45 + 8*M**2*MH**4*N13*N23*N24*N45 - 128*M**2*MH**4*N13*N24*
+ N35*N45 - 64*M**2*MH**4*N13*N24**2*N35 - 64*M**2*MH**4*N13**2*
+ N24*N45 + 8*M**2*MH**4*N14*N23*N24*N35 - 128*M**2*MH**4*N14*N23*
+ N35*N45 - 64*M**2*MH**4*N14*N23**2*N45 - 64*M**2*MH**4*N14**2*
+ N23*N35 + 16*M**2*MH**4*N23*N24*N35*N45 - 18*M**2*N12*N13*N24*
+ N35*X14*X23 - 50*M**2*N12*N13*N24*N35*X14**2 + 50*M**2*N12*N13*
+ N24*X14 - 14*M**2*N12*N13*N24*X23 - 32*M**2*N12*N13*N24**2*N35*
+ X14**2*X23 + 32*M**2*N12*N13*N24**2*X14*X23 - 50*M**2*N12*N13*
+ N35*X14 + 14*M**2*N12*N13*N35*X23 - 18*M**2*N12*N14*N23*N45*X13*
+ X24 - 50*M**2*N12*N14*N23*N45*X13**2 + 50*M**2*N12*N14*N23*X13
+ - 14*M**2*N12*N14*N23*X24 - 32*M**2*N12*N14*N23**2*N45*X13**2*
+ X24 + 32*M**2*N12*N14*N23**2*X13*X24 - 50*M**2*N12*N14*N45*X13
+ + 14*M**2*N12*N14*N45*X24 + 72*M**2*N12*N23*N35*N45*X13*X14 +
+ 36*M**2*N12*N23*N35*N45*X13**2
RES = RES - 36*M**2*N12*N23*N35*X13 - 72*M**2*N12*N23*N35*X14 -
+ 164*M**2*N12*N23*N45*X13 - 14*M**2*N12*N23*N45*X14 - 74*M**2*N12
+ *N23*N45*X24 + 14*M**2*N12*N23*N45**2*X13*X14 + 42*M**2*N12*N23*
+ N45**2*X13*X24 + 46*M**2*N12*N23*N45**2*X13**2 + 118*M**2*N12*
+ N23 - 96*M**2*N12*N23**2*N45*X13*X24 + 32*M**2*N12*N23**2*N45**2
+ *X13**2*X24 + 64*M**2*N12*N23**2*X24 + 72*M**2*N12*N24*N35*N45*
+ X13*X14 + 36*M**2*N12*N24*N35*N45*X14**2 - 14*M**2*N12*N24*N35*
+ X13 - 164*M**2*N12*N24*N35*X14 - 74*M**2*N12*N24*N35*X23 + 14*
+ M**2*N12*N24*N35**2*X13*X14 + 42*M**2*N12*N24*N35**2*X14*X23 +
+ 46*M**2*N12*N24*N35**2*X14**2 - 72*M**2*N12*N24*N45*X13 - 36*
+ M**2*N12*N24*N45*X14 + 118*M**2*N12*N24 - 96*M**2*N12*N24**2*N35
+ *X14*X23 + 32*M**2*N12*N24**2*N35**2*X14**2*X23 + 64*M**2*N12*
+ N24**2*X23 + 108*M**2*N12*N35*N45*X13 + 108*M**2*N12*N35*N45*X14
+ - 118*M**2*N12*N35 + 14*M**2*N12*N35**2*X13 + 46*M**2*N12*
+ N35**2*X14 + 10*M**2*N12*N35**2*X23 - 118*M**2*N12*N45 + 46*M**2
+ *N12*N45**2*X13
RES = RES + 14*M**2*N12*N45**2*X14 + 10*M**2*N12*N45**2*X24 + 4*
+ M**2*N13*N14*N23*N24*X12**2 + 8*M**2*N13*N14*N23*X12 + 4*M**2*
+ N13*N14*N23*X24 + 8*M**2*N13*N14*N24*X12 + 4*M**2*N13*N14*N24*
+ X23 + 4*M**2*N13*N14*N35*N45*X12**2 - 4*M**2*N13*N14*N35*X12 + 4
+ *M**2*N13*N14*N35*X23 - 4*M**2*N13*N14*N45*X12 + 4*M**2*N13*N14*
+ N45*X24 + 16*M**2*N13*N14 + 8*M**2*N13*N23*N24*X12 + 4*M**2*N13*
+ N23*N24*X14 + 8*M**2*N13*N23*N45*X12 + 8*M**2*N13*N23*N45*X14 +
+ 8*M**2*N13*N23*N45*X24 - 4*M**2*N13*N23*N45**2*X12*X14 - 4*M**2*
+ N13*N23*N45**2*X12*X24 - 4*M**2*N13*N23*N45**2*X12**2 + 8*M**2*
+ N13*N23 - 100*M**2*N13*N24*N35*X12 - 118*M**2*N13*N24*N35*X14 -
+ 50*M**2*N13*N24*N35*X23 - 32*M**2*N13*N24*N45*X12 - 32*M**2*N13*
+ N24*N45*X14 + 72*M**2*N13*N24 - 64*M**2*N13*N24**2*N35*X12*X14
+ - 32*M**2*N13*N24**2*N35*X12*X23 - 32*M**2*N13*N24**2*N35*
+ X12**2 - 64*M**2*N13*N24**2*N35*X14*X23 - 32*M**2*N13*N24**2*N35
+ *X14**2 + 32*M**2*N13*N24**2*X23 + 40*M**2*N13*N35*N45*X12 + 36*
+ M**2*N13*N35*N45*X14
RES = RES - 108*M**2*N13*N35 + 28*M**2*N13*N45 - 4*M**2*N13*
+ N45**2*X12 - 4*M**2*N13*N45**2*X14 - 4*M**2*N13*N45**2*X24 - 32*
+ M**2*N13**2*N24*X12 - 32*M**2*N13**2*N24*X23 - 32*M**2*N13**2*
+ N45*X24 - 64*M**2*N13**2 + 8*M**2*N14*N23*N24*X12 + 4*M**2*N14*
+ N23*N24*X13 - 32*M**2*N14*N23*N35*X12 - 32*M**2*N14*N23*N35*X13
+ - 100*M**2*N14*N23*N45*X12 - 118*M**2*N14*N23*N45*X13 - 50*M**2
+ *N14*N23*N45*X24 + 72*M**2*N14*N23 - 64*M**2*N14*N23**2*N45*X12*
+ X13 - 32*M**2*N14*N23**2*N45*X12*X24 - 32*M**2*N14*N23**2*N45*
+ X12**2 - 64*M**2*N14*N23**2*N45*X13*X24 - 32*M**2*N14*N23**2*N45
+ *X13**2 + 32*M**2*N14*N23**2*X24 + 8*M**2*N14*N24*N35*X12 + 8*
+ M**2*N14*N24*N35*X13 + 8*M**2*N14*N24*N35*X23 - 4*M**2*N14*N24*
+ N35**2*X12*X13 - 4*M**2*N14*N24*N35**2*X12*X23 - 4*M**2*N14*N24*
+ N35**2*X12**2 + 8*M**2*N14*N24 + 40*M**2*N14*N35*N45*X12 + 36*
+ M**2*N14*N35*N45*X13 + 28*M**2*N14*N35 - 4*M**2*N14*N35**2*X12
+ - 4*M**2*N14*N35**2*X13 - 4*M**2*N14*N35**2*X23 - 108*M**2*N14*
+ N45
RES = RES - 32*M**2*N14**2*N23*X12 - 32*M**2*N14**2*N23*X24 - 32*
+ M**2*N14**2*N35*X23 - 64*M**2*N14**2 + 16*M**2*N23*N24*N35*N45*
+ X12*X13 + 16*M**2*N23*N24*N35*N45*X12*X14 + 12*M**2*N23*N24*N35*
+ N45*X12**2 + 16*M**2*N23*N24*N35*N45*X13*X14 + 4*M**2*N23*N24*
+ N35*N45*X13**2 + 4*M**2*N23*N24*N35*N45*X14**2 - 8*M**2*N23*N24*
+ N35*X12 - 8*M**2*N23*N24*N35*X14 - 8*M**2*N23*N24*N45*X12 - 8*
+ M**2*N23*N24*N45*X13 + 16*M**2*N23*N24 + 44*M**2*N23*N35*N45*X12
+ + 80*M**2*N23*N35*N45*X13 + 40*M**2*N23*N35*N45*X14 - 108*M**2*
+ N23*N35 - 206*M**2*N23*N45 + 84*M**2*N23*N45**2*X12 + 134*M**2*
+ N23*N45**2*X13 + 10*M**2*N23*N45**2*X14 + 38*M**2*N23*N45**2*X24
+ - 128*M**2*N23**2*N45*X12 - 128*M**2*N23**2*N45*X13 - 32*M**2*
+ N23**2*N45*X14 - 96*M**2*N23**2*N45*X24 + 64*M**2*N23**2*N45**2*
+ X12*X13 + 32*M**2*N23**2*N45**2*X12*X24 + 32*M**2*N23**2*N45**2*
+ X12**2 + 64*M**2*N23**2*N45**2*X13*X24 + 32*M**2*N23**2*N45**2*
+ X13**2 + 32*M**2*N23**2 + 44*M**2*N24*N35*N45*X12 + 40*M**2*N24*
+ N35*N45*X13
RES = RES + 80*M**2*N24*N35*N45*X14 - 206*M**2*N24*N35 + 84*M**2*
+ N24*N35**2*X12 + 10*M**2*N24*N35**2*X13 + 134*M**2*N24*N35**2*
+ X14 + 38*M**2*N24*N35**2*X23 - 108*M**2*N24*N45 - 128*M**2*
+ N24**2*N35*X12 - 32*M**2*N24**2*N35*X13 - 128*M**2*N24**2*N35*
+ X14 - 96*M**2*N24**2*N35*X23 + 64*M**2*N24**2*N35**2*X12*X14 +
+ 32*M**2*N24**2*N35**2*X12*X23 + 32*M**2*N24**2*N35**2*X12**2 +
+ 64*M**2*N24**2*N35**2*X14*X23 + 32*M**2*N24**2*N35**2*X14**2 +
+ 32*M**2*N24**2 + 216*M**2*N35*N45 + 56*M**2*N35**2 + 56*M**2*
+ N45**2 - 16*M**4*MH**2*N13*N14*N23*N24 - 16*M**4*MH**2*N13*N14*
+ N23*N45 - 16*M**4*MH**2*N13*N14*N24*N35 - 16*M**4*MH**2*N13*N14*
+ N35*N45 - 16*M**4*MH**2*N13*N23*N24*N45 - 16*M**4*MH**2*N13*N23*
+ N45**2 + 128*M**4*MH**2*N13*N24*N35*N45 + 128*M**4*MH**2*N13*
+ N24**2*N35 + 128*M**4*MH**2*N13**2*N24*N45 + 64*M**4*MH**2*
+ N13**2*N24**2 + 64*M**4*MH**2*N13**2*N45**2 - 16*M**4*MH**2*N14*
+ N23*N24*N35 + 128*M**4*MH**2*N14*N23*N35*N45 + 128*M**4*MH**2*
+ N14*N23**2*N45
RES = RES - 16*M**4*MH**2*N14*N24*N35**2 + 128*M**4*MH**2*N14**2*
+ N23*N35 + 64*M**4*MH**2*N14**2*N23**2 + 64*M**4*MH**2*N14**2*
+ N35**2 - 16*M**4*MH**2*N23*N24*N35*N45 + 64*M**4*MH**2*N23**2*
+ N45**2 + 64*M**4*MH**2*N24**2*N35**2 - MH**2*N12*N13*N24*N35*X14
+ *X23 - 5*MH**2*N12*N13*N24*N35*X14**2 - 31*MH**2*N12*N13*N24*X14
+ - 49*MH**2*N12*N13*N24*X23 + 31*MH**2*N12*N13*N35*X14 - 23*
+ MH**2*N12*N13*N35*X23 - 72*MH**2*N12*N13 - MH**2*N12*N14*N23*N45
+ *X13*X24 - 5*MH**2*N12*N14*N23*N45*X13**2 - 31*MH**2*N12*N14*N23
+ *X13 - 49*MH**2*N12*N14*N23*X24 + 31*MH**2*N12*N14*N45*X13 - 23*
+ MH**2*N12*N14*N45*X24 - 72*MH**2*N12*N14 - 18*MH**2*N12*N23*N35*
+ N45*X13*X14 - 40*MH**2*N12*N23*N35*N45*X13**2 + 4*MH**2*N12*N23*
+ N35*X13 + 18*MH**2*N12*N23*N35*X14 - 41*MH**2*N12*N23*N45*X13 +
+ 23*MH**2*N12*N23*N45**2*X13**2 - 54*MH**2*N12*N23 - 18*MH**2*N12
+ *N24*N35*N45*X13*X14 - 40*MH**2*N12*N24*N35*N45*X14**2 - 41*
+ MH**2*N12*N24*N35*X14 + 23*MH**2*N12*N24*N35**2*X14**2 + 18*
+ MH**2*N12*N24*N45*X13
RES = RES + 4*MH**2*N12*N24*N45*X14 - 54*MH**2*N12*N24 - 166*
+ MH**2*N12*N35*N45*X13 - 166*MH**2*N12*N35*N45*X14 + 162*MH**2*
+ N12*N35 + 95*MH**2*N12*N35**2*X14 + 162*MH**2*N12*N45 + 95*MH**2
+ *N12*N45**2*X13 - 14*MH**2*N12**2*N13*N24*N35*X14**2*X23 + 14*
+ MH**2*N12**2*N13*N24*X14*X23 - 14*MH**2*N12**2*N13*N35*X14*X23
+ - 14*MH**2*N12**2*N14*N23*N45*X13**2*X24 + 14*MH**2*N12**2*N14*
+ N23*X13*X24 - 14*MH**2*N12**2*N14*N45*X13*X24 + 36*MH**2*N12**2*
+ N23*N35*N45*X13**2*X14 - 36*MH**2*N12**2*N23*N35*X13*X14 - 36*
+ MH**2*N12**2*N23*N45*X13**2 + 36*MH**2*N12**2*N23*X13 + 36*MH**2
+ *N12**2*N24*N35*N45*X13*X14**2 - 36*MH**2*N12**2*N24*N35*X14**2
+ - 36*MH**2*N12**2*N24*N45*X13*X14 + 36*MH**2*N12**2*N24*X14 -
+ 72*MH**2*N12**2*N35*N45*X13*X14 - 36*MH**2*N12**2*N35*X14 + 72*
+ MH**2*N12**2*N35**2*X14**2 - 36*MH**2*N12**2*N45*X13 + 72*MH**2*
+ N12**2*N45**2*X13**2 - 6*MH**2*N13*N14*N23*N24*X12**2 + 2*MH**2*
+ N13*N14*N23*N45*X12*X24 + 2*MH**2*N13*N14*N23*N45*X12**2 - 14*
+ MH**2*N13*N14*N23*X12
RES = RES - 8*MH**2*N13*N14*N23*X24 + 2*MH**2*N13*N14*N24*N35*X12
+ *X23 + 2*MH**2*N13*N14*N24*N35*X12**2 - 14*MH**2*N13*N14*N24*X12
+ - 8*MH**2*N13*N14*N24*X23 - 6*MH**2*N13*N14*N35*N45*X12**2 + 8*
+ MH**2*N13*N14*N35*X12 - 4*MH**2*N13*N14*N35*X23 + 8*MH**2*N13*
+ N14*N45*X12 - 4*MH**2*N13*N14*N45*X24 - 24*MH**2*N13*N14 + 2*
+ MH**2*N13*N23*N24*N45*X12*X14 + 2*MH**2*N13*N23*N24*N45*X12**2
+ - 14*MH**2*N13*N23*N24*X12 - 8*MH**2*N13*N23*N24*X14 + 8*MH**2*
+ N13*N23*N45*X12 - 2*MH**2*N13*N23*N45**2*X12**2 - 18*MH**2*N13*
+ N23 + 47*MH**2*N13*N24*N35*X12 + 24*MH**2*N13*N24*N35*X14 + 31*
+ MH**2*N13*N24*N35*X23 + 18*MH**2*N13*N24*N45*X12 + 18*MH**2*N13*
+ N24*N45*X14 - 77*MH**2*N13*N24 - 60*MH**2*N13*N35*N45*X12 - 54*
+ MH**2*N13*N35*N45*X14 + 105*MH**2*N13*N35 + 60*MH**2*N13*N45 +
+ 14*MH**2*N13*N45**2*X12 + 2*MH**2*N14*N23*N24*N35*X12*X13 + 2*
+ MH**2*N14*N23*N24*N35*X12**2 - 14*MH**2*N14*N23*N24*X12 - 8*
+ MH**2*N14*N23*N24*X13 + 18*MH**2*N14*N23*N35*X12 + 18*MH**2*N14*
+ N23*N35*X13
RES = RES + 47*MH**2*N14*N23*N45*X12 + 24*MH**2*N14*N23*N45*X13
+ + 31*MH**2*N14*N23*N45*X24 - 77*MH**2*N14*N23 + 8*MH**2*N14*N24
+ *N35*X12 - 2*MH**2*N14*N24*N35**2*X12**2 - 18*MH**2*N14*N24 - 60
+ *MH**2*N14*N35*N45*X12 - 54*MH**2*N14*N35*N45*X13 + 60*MH**2*N14
+ *N35 + 14*MH**2*N14*N35**2*X12 + 105*MH**2*N14*N45 - 22*MH**2*
+ N23*N24*N35*N45*X12*X13 - 22*MH**2*N23*N24*N35*N45*X12*X14 - 14*
+ MH**2*N23*N24*N35*N45*X12**2 - 20*MH**2*N23*N24*N35*N45*X13*X14
+ - 8*MH**2*N23*N24*N35*N45*X13**2 - 8*MH**2*N23*N24*N35*N45*
+ X14**2 + 10*MH**2*N23*N24*N35*X12 + 4*MH**2*N23*N24*N35*X13 + 8*
+ MH**2*N23*N24*N35*X14 + 10*MH**2*N23*N24*N45*X12 + 8*MH**2*N23*
+ N24*N45*X13 + 4*MH**2*N23*N24*N45*X14 - 20*MH**2*N23*N24 - 14*
+ MH**2*N23*N35*N45*X12 - 54*MH**2*N23*N35*N45*X13 - 8*MH**2*N23*
+ N35*N45*X14 + 42*MH**2*N23*N35 + 35*MH**2*N23*N45 + 3*MH**2*N23*
+ N45**2*X12 + 12*MH**2*N23*N45**2*X13 - 14*MH**2*N24*N35*N45*X12
+ - 8*MH**2*N24*N35*N45*X13 - 54*MH**2*N24*N35*N45*X14 + 35*MH**2
+ *N24*N35
RES = RES + 3*MH**2*N24*N35**2*X12 + 12*MH**2*N24*N35**2*X14 + 42
+ *MH**2*N24*N45 - 232*MH**2*N35*N45 + 57*MH**2*N35**2 + 57*MH**2*
+ N45**2 - 18*MH**4*N12*N13*N24*N35*X14 + 54*MH**4*N12*N13*N24 +
+ 18*MH**4*N12*N13*N35 - 18*MH**4*N12*N14*N23*N45*X13 + 54*MH**4*
+ N12*N14*N23 + 18*MH**4*N12*N14*N45 + 72*MH**4*N12*N23*N35*N45*
+ X13 - 36*MH**4*N12*N23*N35 + 72*MH**4*N12*N24*N35*N45*X14 - 36*
+ MH**4*N12*N24*N45 - 144*MH**4*N12*N35*N45 + 8*MH**4*N13*N14*N23*
+ N24*X12 - 2*MH**4*N13*N14*N23*N45*X12 + 10*MH**4*N13*N14*N23 - 2
+ *MH**4*N13*N14*N24*N35*X12 + 10*MH**4*N13*N14*N24 - 8*MH**4*N13*
+ N14*N35*N45*X12 + 6*MH**4*N13*N14*N35 + 6*MH**4*N13*N14*N45 - 2*
+ MH**4*N13*N23*N24*N45*X12 + 10*MH**4*N13*N23*N24 - 16*MH**4*N13*
+ N24*N35*N45*X12 - 16*MH**4*N13*N24*N35*N45*X14 - 52*MH**4*N13*
+ N24*N35 - 18*MH**4*N13*N24*N45 - 88*MH**4*N13*N35*N45 - 2*MH**4*
+ N14*N23*N24*N35*X12 + 10*MH**4*N14*N23*N24 - 16*MH**4*N14*N23*
+ N35*N45*X12 - 16*MH**4*N14*N23*N35*N45*X13 - 18*MH**4*N14*N23*
+ N35
RES = RES - 52*MH**4*N14*N23*N45 - 88*MH**4*N14*N35*N45 + 12*
+ MH**4*N23*N24*N35*N45*X12 + 10*MH**4*N23*N24*N35*N45*X13 + 10*
+ MH**4*N23*N24*N35*N45*X14 - 4*MH**4*N23*N24*N35 - 4*MH**4*N23*
+ N24*N45 + 10*MH**4*N23*N35*N45 + 10*MH**4*N24*N35*N45 - 4*MH**6*
+ N13*N14*N23*N24 - 4*MH**6*N13*N14*N35*N45 + 32*MH**6*N13*N24*N35
+ *N45 + 32*MH**6*N14*N23*N35*N45 - 4*MH**6*N23*N24*N35*N45 - 2*
+ N12*N13*N24*N35*X14*X23**2 + 28*N12*N13*N24*N35*X14**2*X23 + 7*
+ N12*N13*N24*N35*X14**3 + 15*N12*N13*N24*X14*X23 + 11*N12*N13*N24
+ *X14**2 + 27*N12*N13*N24*X23**2 - 15*N12*N13*N35*X14*X23 + 25*
+ N12*N13*N35*X14**2 + 9*N12*N13*N35*X23**2 + 18*N12*N13*X14 + 72*
+ N12*N13*X23 + 36*N12*N13*X24 - 2*N12*N14*N23*N45*X13*X24**2 + 28
+ *N12*N14*N23*N45*X13**2*X24 + 7*N12*N14*N23*N45*X13**3 + 15*N12*
+ N14*N23*X13*X24 + 11*N12*N14*N23*X13**2 + 27*N12*N14*N23*X24**2
+ - 15*N12*N14*N45*X13*X24 + 25*N12*N14*N45*X13**2 + 9*N12*N14*
+ N45*X24**2 + 18*N12*N14*X13 + 36*N12*N14*X23 + 72*N12*N14*X24 -
+ 16*N12*N23*N35*N45*X13**2*X14
RES = RES + 2*N12*N23*N35*N45*X13**3 + 34*N12*N23*N35*X13*X14 +
+ 16*N12*N23*N35*X13**2 + 9*N12*N23*N45*X13*X14 + 45*N12*N23*N45*
+ X13*X24 + 98*N12*N23*N45*X13**2 + 9*N12*N23*N45**2*X13**2*X14 -
+ 35*N12*N23*N45**2*X13**2*X24 - 23*N12*N23*N45**2*X13**3 - 21*N12
+ *N23*X13 + 18*N12*N23*X14 + 33*N12*N23*X24 - 16*N12*N24*N35*N45*
+ X13*X14**2 + 2*N12*N24*N35*N45*X14**3 + 9*N12*N24*N35*X13*X14 +
+ 45*N12*N24*N35*X14*X23 + 98*N12*N24*N35*X14**2 + 9*N12*N24*
+ N35**2*X13*X14**2 - 35*N12*N24*N35**2*X14**2*X23 - 23*N12*N24*
+ N35**2*X14**3 + 34*N12*N24*N45*X13*X14 + 16*N12*N24*N45*X14**2
+ + 18*N12*N24*X13 - 21*N12*N24*X14 + 33*N12*N24*X23 - 140*N12*
+ N35*N45*X13*X14 - 52*N12*N35*N45*X13**2 - 52*N12*N35*N45*X14**2
+ + 36*N12*N35*X13 + 129*N12*N35*X14 - 33*N12*N35*X23 + 9*N12*
+ N35**2*X13*X14 - 12*N12*N35**2*X14*X23 - 23*N12*N35**2*X14**2 +
+ 129*N12*N45*X13 + 36*N12*N45*X14 - 33*N12*N45*X24 + 9*N12*N45**2
+ *X13*X14 - 12*N12*N45**2*X13*X24 - 23*N12*N45**2*X13**2 + 7*
+ N12**2*N13*N24*N35*X14**2*X23**2
RES = RES + 7*N12**2*N13*N24*N35*X14**3*X23 - 7*N12**2*N13*N24*
+ X14*X23**2 - 7*N12**2*N13*N24*X14**2*X23 + 7*N12**2*N13*N35*X14*
+ X23**2 + 7*N12**2*N13*N35*X14**2*X23 + 7*N12**2*N14*N23*N45*
+ X13**2*X24**2 + 7*N12**2*N14*N23*N45*X13**3*X24 - 7*N12**2*N14*
+ N23*X13*X24**2 - 7*N12**2*N14*N23*X13**2*X24 + 7*N12**2*N14*N45*
+ X13*X24**2 + 7*N12**2*N14*N45*X13**2*X24 - 18*N12**2*N23*N35*N45
+ *X13**2*X14**2 - 18*N12**2*N23*N35*N45*X13**3*X14 + 18*N12**2*
+ N23*N35*X13*X14**2 + 18*N12**2*N23*N35*X13**2*X14 + 18*N12**2*
+ N23*N45*X13**2*X14 + 53*N12**2*N23*N45*X13**2*X24 + 18*N12**2*
+ N23*N45*X13**3 - 23*N12**2*N23*N45**2*X13**3*X24 - 18*N12**2*N23
+ *X13*X14 - 30*N12**2*N23*X13*X24 - 18*N12**2*N23*X13**2 - 18*
+ N12**2*N24*N35*N45*X13*X14**3 - 18*N12**2*N24*N35*N45*X13**2*
+ X14**2 + 18*N12**2*N24*N35*X13*X14**2 + 53*N12**2*N24*N35*X14**2
+ *X23 + 18*N12**2*N24*N35*X14**3 - 23*N12**2*N24*N35**2*X14**3*
+ X23 + 18*N12**2*N24*N45*X13*X14**2 + 18*N12**2*N24*N45*X13**2*
+ X14
RES = RES - 18*N12**2*N24*X13*X14 - 30*N12**2*N24*X14*X23 - 18*
+ N12**2*N24*X14**2 - 36*N12**2*N35*N45*X13*X14**2 - 36*N12**2*N35
+ *N45*X13**2*X14 + 18*N12**2*N35*X13*X14 + 30*N12**2*N35*X14*X23
+ + 18*N12**2*N35*X14**2 - 23*N12**2*N35**2*X14**2*X23 + 18*
+ N12**2*N45*X13*X14 + 30*N12**2*N45*X13*X24 + 18*N12**2*N45*
+ X13**2 - 23*N12**2*N45**2*X13**2*X24 + 2*N13*N14*N23*N24*X12**3
+ + 6*N13*N14*N23*X12*X24 + 6*N13*N14*N23*X12**2 + 2*N13*N14*N23*
+ X24**2 + 6*N13*N14*N24*X12*X23 + 6*N13*N14*N24*X12**2 + 2*N13*
+ N14*N24*X23**2 - 2*N13*N14*N35*N45*X12**3 - 2*N13*N14*N35*X12*
+ X23 + 2*N13*N14*N35*X12**2 + 2*N13*N14*N35*X23**2 - 2*N13*N14*
+ N45*X12*X24 + 2*N13*N14*N45*X12**2 + 2*N13*N14*N45*X24**2 + 12*
+ N13*N14*X12 + 12*N13*N14*X23 + 12*N13*N14*X24 + 6*N13*N23*N24*
+ X12*X14 + 6*N13*N23*N24*X12**2 + 2*N13*N23*N24*X14**2 - 4*N13*
+ N23*N45*X12*X14 - 4*N13*N23*N45*X12*X24 - 4*N13*N23*N45*X12**2
+ + 2*N13*N23*N45**2*X12**2*X14 + 2*N13*N23*N45**2*X12**2*X24 + 2
+ *N13*N23*N45**2*X12**3
RES = RES + 14*N13*N23*X12 + 8*N13*N23*X14 + 8*N13*N23*X24 + 5*
+ N13*N24*N35*X12*X14 - 18*N13*N24*N35*X12*X23 - 9*N13*N24*N35*
+ X12**2 + 3*N13*N24*N35*X14*X23 + 21*N13*N24*N35*X14**2 - 9*N13*
+ N24*N35*X23**2 + 39*N13*N24*X12 + 28*N13*N24*X14 + 60*N13*N24*
+ X23 - 38*N13*N35*N45*X12*X14 - 22*N13*N35*N45*X12**2 - 18*N13*
+ N35*N45*X14**2 + 13*N13*N35*X12 + 52*N13*N35*X14 - 40*N13*N35*
+ X23 + 18*N13*N45*X12 + 16*N13*N45*X14 - 20*N13*N45*X24 + 2*N13*
+ N45**2*X12*X14 + 2*N13*N45**2*X12*X24 + 2*N13*N45**2*X12**2 + 28
+ *N13 + 6*N14*N23*N24*X12*X13 + 6*N14*N23*N24*X12**2 + 2*N14*N23*
+ N24*X13**2 + 5*N14*N23*N45*X12*X13 - 18*N14*N23*N45*X12*X24 - 9*
+ N14*N23*N45*X12**2 + 3*N14*N23*N45*X13*X24 + 21*N14*N23*N45*
+ X13**2 - 9*N14*N23*N45*X24**2 + 39*N14*N23*X12 + 28*N14*N23*X13
+ + 60*N14*N23*X24 - 4*N14*N24*N35*X12*X13 - 4*N14*N24*N35*X12*
+ X23 - 4*N14*N24*N35*X12**2 + 2*N14*N24*N35**2*X12**2*X13 + 2*N14
+ *N24*N35**2*X12**2*X23 + 2*N14*N24*N35**2*X12**3 + 14*N14*N24*
+ X12
RES = RES + 8*N14*N24*X13 + 8*N14*N24*X23 - 38*N14*N35*N45*X12*
+ X13 - 22*N14*N35*N45*X12**2 - 18*N14*N35*N45*X13**2 + 18*N14*N35
+ *X12 + 16*N14*N35*X13 - 20*N14*N35*X23 + 2*N14*N35**2*X12*X13 +
+ 2*N14*N35**2*X12*X23 + 2*N14*N35**2*X12**2 + 13*N14*N45*X12 + 52
+ *N14*N45*X13 - 40*N14*N45*X24 + 28*N14 + 24*N23*N24*N35*N45*X12*
+ X13*X14 + 10*N23*N24*N35*N45*X12*X13**2 + 10*N23*N24*N35*N45*X12
+ *X14**2 + 14*N23*N24*N35*N45*X12**2*X13 + 14*N23*N24*N35*N45*
+ X12**2*X14 + 6*N23*N24*N35*N45*X12**3 + 10*N23*N24*N35*N45*X13*
+ X14**2 + 10*N23*N24*N35*N45*X13**2*X14 + 2*N23*N24*N35*N45*
+ X13**3 + 2*N23*N24*N35*N45*X14**3 - 4*N23*N24*N35*X12*X13 - 8*
+ N23*N24*N35*X12*X14 - 4*N23*N24*N35*X12**2 - 4*N23*N24*N35*X13*
+ X14 - 4*N23*N24*N35*X14**2 - 8*N23*N24*N45*X12*X13 - 4*N23*N24*
+ N45*X12*X14 - 4*N23*N24*N45*X12**2 - 4*N23*N24*N45*X13*X14 - 4*
+ N23*N24*N45*X13**2 + 16*N23*N24*X12 + 10*N23*N24*X13 + 10*N23*
+ N24*X14 + 14*N23*N35*N45*X12*X13 + 8*N23*N35*N45*X12*X14 + 6*N23
+ *N35*N45*X12**2
RES = RES + 10*N23*N35*N45*X13*X14 + 10*N23*N35*N45*X13**2 + 2*
+ N23*N35*N45*X14**2 - 4*N23*N35*X12 + 12*N23*N35*X13 - 4*N23*N35*
+ X14 + 6*N23*N45*X12 + 74*N23*N45*X13 - 13*N23*N45*X14 + 4*N23*
+ N45*X24 - 15*N23*N45**2*X12*X13 + 13*N23*N45**2*X12*X14 - N23*
+ N45**2*X12*X24 - N23*N45**2*X12**2 + 20*N23*N45**2*X13*X14 - 15*
+ N23*N45**2*X13*X24 - 35*N23*N45**2*X13**2 + 43*N23 + 8*N24*N35*
+ N45*X12*X13 + 14*N24*N35*N45*X12*X14 + 6*N24*N35*N45*X12**2 + 10
+ *N24*N35*N45*X13*X14 + 2*N24*N35*N45*X13**2 + 10*N24*N35*N45*
+ X14**2 + 6*N24*N35*X12 - 13*N24*N35*X13 + 74*N24*N35*X14 + 4*N24
+ *N35*X23 + 13*N24*N35**2*X12*X13 - 15*N24*N35**2*X12*X14 - N24*
+ N35**2*X12*X23 - N24*N35**2*X12**2 + 20*N24*N35**2*X13*X14 - 15*
+ N24*N35**2*X14*X23 - 35*N24*N35**2*X14**2 - 4*N24*N45*X12 - 4*
+ N24*N45*X13 + 12*N24*N45*X14 + 43*N24 - 108*N35*N45*X12 - 138*
+ N35*N45*X13 - 138*N35*N45*X14 + 95*N35 - 3*N35**2*X12 + 11*
+ N35**2*X13 - 12*N35**2*X14 - 3*N35**2*X23 + 95*N45 - 3*N45**2*
+ X12
RES = RES - 12*N45**2*X13 + 11*N45**2*X14 - 3*N45**2*X24
DAGG = RES
RETURN
END
DOUBLE PRECISION FUNCTION DLUMQQ(TAU,X,Q)
C--Q QBAR-LUMINOSITY
IMPLICIT DOUBLE PRECISION (A-H,O-Z)
DIMENSION PDF1(-6:6),PDF2(-6:6)
COMMON/COLLIDER/ICOLL
CALL STRUC(X,Q,PDF1)
CALL STRUC(TAU/X,Q,PDF2)
DL = 0
DO I=1,5
DL = DL + PDF1(I)*PDF2(-ICOLL*I) + PDF1(-I)*PDF2(ICOLL*I)
ENDDO
DLUMQQ = DL/TAU
RETURN
END
DOUBLE PRECISION FUNCTION DLUMGG(TAU,X,Q)
C--G G-LUMINOSITY
IMPLICIT DOUBLE PRECISION (A-H,O-Z)
DIMENSION PDF1(-6:6),PDF2(-6:6)
CALL STRUC(X,Q,PDF1)
CALL STRUC(TAU/X,Q,PDF2)
DLUMGG = PDF1(0)*PDF2(0)/TAU
RETURN
END
DOUBLE PRECISION FUNCTION ALPHAS(Q,N)
C--ALPHA_S: Q = SCALE, N = 1 -> LO, N = 2 -> NLO
IMPLICIT DOUBLE PRECISION (A-H,O-Z)
DIMENSION XLB(6)
COMMON/ALSLAM/XLB1(6),XLB2(6)
COMMON/ALS/XLAMBDA,AMC,AMB,AMT,N0
B0(NF)=33.D0-2.D0*NF
B1(NF)=6.D0*(153.D0-19.D0*NF)/B0(NF)**2
ALS1(NF,X)=12.D0*PI/(B0(NF)*DLOG(X**2/XLB(NF)**2))
ALS2(NF,X)=12.D0*PI/(B0(NF)*DLOG(X**2/XLB(NF)**2))
. *(1.D0-B1(NF)*DLOG(DLOG(X**2/XLB(NF)**2))
. /DLOG(X**2/XLB(NF)**2))
PI=4.D0*DATAN(1.D0)
IF(N.EQ.1)THEN
DO 1 I=1,6
XLB(I)=XLB1(I)
1 CONTINUE
ELSE
DO 2 I=1,6
XLB(I)=XLB2(I)
2 CONTINUE
ENDIF
IF(Q.LT.AMC)THEN
NF=3
ELSEIF(Q.LE.AMB)THEN
NF=4
ELSEIF(Q.LE.AMT)THEN
NF=5
ELSE
NF=6
ENDIF
IF(N.EQ.1)THEN
ALPHAS=ALS1(NF,Q)
ELSE
ALPHAS=ALS2(NF,Q)
ENDIF
RETURN
END
SUBROUTINE ALSINI(ACC)
C--ALPHA_S: INITIALIZATION OF LAMBDA_NF, ACC = ACCURAY AT THRESHOLDS
IMPLICIT DOUBLE PRECISION (A-H,O-Z)
DIMENSION XLB(6)
COMMON/ALSLAM/XLB1(6),XLB2(6)
COMMON/ALS/XLAMBDA,AMC,AMB,AMT,N0
B0(NF)=33.D0-2.D0*NF
B1(NF)=6.D0*(153.D0-19.D0*NF)/B0(NF)**2
ALS2(NF,X)=12.D0*PI/(B0(NF)*DLOG(X**2/XLB(NF)**2))
. *(1.D0-B1(NF)*DLOG(DLOG(X**2/XLB(NF)**2))
. /DLOG(X**2/XLB(NF)**2))
PI=4.D0*DATAN(1.D0)
XLB1(1)=0D0
XLB1(2)=0D0
XLB2(1)=0D0
XLB2(2)=0D0
IF(N0.EQ.4)THEN
XLB(4)=XLAMBDA
XLB(5)=XLB(4)*(XLB(4)/AMB)**(2.D0/23.D0)
ELSEIF(N0.EQ.5)THEN
XLB(5)=XLAMBDA
XLB(4)=XLB(5)*(XLB(5)/AMB)**(-2.D0/25.D0)
ENDIF
XLB(3)=XLB(4)*(XLB(4)/AMC)**(-2.D0/27.D0)
XLB(6)=XLB(5)*(XLB(5)/AMT)**(2.D0/21.D0)
DO 1 I=1,6
XLB1(I)=XLB(I)
1 CONTINUE
IF(N0.EQ.4)THEN
XLB(4)=XLAMBDA
XLB(5)=XLB(4)*(XLB(4)/AMB)**(2.D0/23.D0)
. *(2.D0*DLOG(AMB/XLB(4)))**(-963.D0/13225.D0)
XLB(5)=XITER(AMB,XLB(4),4,XLB(5),5,ACC)
ELSEIF(N0.EQ.5)THEN
XLB(5)=XLAMBDA
XLB(4)=XLB(5)*(XLB(5)/AMB)**(-2.D0/25.D0)
. *(2.D0*DLOG(AMB/XLB(5)))**(963.D0/14375.D0)
XLB(4)=XITER(AMB,XLB(5),5,XLB(4),4,ACC)
ENDIF
XLB(3)=XLB(4)*(XLB(4)/AMC)**(-2.D0/27.D0)
. *(2.D0*DLOG(AMC/XLB(4)))**(107.D0/2025.D0)
XLB(3)=XITER(AMC,XLB(4),4,XLB(3),3,ACC)
XLB(6)=XLB(5)*(XLB(5)/AMT)**(2.D0/21.D0)
. *(2.D0*DLOG(AMT/XLB(5)))**(-321.D0/3381.D0)
XLB(6)=XITER(AMT,XLB(5),5,XLB(6),6,ACC)
DO 2 I=1,6
XLB2(I)=XLB(I)
2 CONTINUE
RETURN
END
DOUBLE PRECISION FUNCTION XITER(Q,XLB1,NF1,XLB,NF2,ACC)
C--ALPHA_S: ITERATION FOR ALPHA_S(M_Z) -> LAMBDA_5
IMPLICIT DOUBLE PRECISION (A-H,O-Z)
B0(NF)=33.D0-2.D0*NF
B1(NF)=6.D0*(153.D0-19.D0*NF)/B0(NF)**2
ALS2(NF,X,XLB)=12.D0*PI/(B0(NF)*DLOG(X**2/XLB**2))
. *(1.D0-B1(NF)*DLOG(DLOG(X**2/XLB**2))
. /DLOG(X**2/XLB**2))
AA(NF)=12D0*PI/B0(NF)
BB(NF)=B1(NF)/AA(NF)
XIT(A,B,X)=A/2.D0*(1D0+DSQRT(1D0-4D0*B*DLOG(X)))
PI=4.D0*DATAN(1.D0)
XLB2=XLB
IF(ACC.GE.1.D0) GOTO 111
II=0
1 II=II+1
X=DLOG(Q**2/XLB2**2)
ALP=ALS2(NF1,Q,XLB1)
A=AA(NF2)/ALP
B=BB(NF2)*ALP
XX=XIT(A,B,X)
XLB2=Q*DEXP(-XX/2.D0)
Y1=ALS2(NF1,Q,XLB1)
Y2=ALS2(NF2,Q,XLB2)
DY=DABS(Y2-Y1)/Y1
IF(DY.GE.ACC) GOTO 1
111 XITER=XLB2
RETURN
END
SUBROUTINE VEGASN(FXN,ACC,NDIM,NPOINT,NITT,NPRN,IGRAPH)
C--VEGAS: WARMUP AND MAIN INTEGRATION
IMPLICIT DOUBLE PRECISION (A-H,O-Z)
EXTERNAL FXN
NPOINT0=NPOINT/10
NITT0=2*NITT
CALL VEGAS(FXN,ACC,NDIM,NPOINT0,NITT0,NPRN,IGRAPH)
CALL VEGAS1(FXN,ACC,NDIM,NPOINT,NITT,NPRN,IGRAPH)
RETURN
END
C%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
SUBROUTINE VEGAS(FXN,BCC,NDIM,NCALL,ITMX,NPRN,IGRAPH)
C--VEGAS SUBROUTINE
IMPLICIT DOUBLE PRECISION (A-H,O-Z)
COMMON/RESU/RES
COMMON/VEGOUT/NV
COMMON/BVEG2/NDO,IT,SI,SI2,SWGT,SCHI,XI(50,10),SCALLS
1,D(50,10),DI(50,10),NXI(50,10)
DIMENSION XIN(50),R(50),DX(10),IA(10),KG(10),DT(10)
DIMENSION XL(10),XU(10),QRAN(10),X(10)
COMMON/RESULT/S1,S2,S3,S4
EXTERNAL FXN
DATA XL,XU/10*0.D0,10*1.D0/
DATA NDMX/50/,ALPH/1.5D0/,ONE/1.D0/,MDS/1/
IPR=1
IF(NPRN.GT.0)IPR=0
NDO=1
DO 1 J=1,NDIM
1 XI(1,J)=ONE
ENTRY VEGAS1(FXN,BCC,NDIM,NCALL,ITMX,NPRN,IGRAPH)
NOW=IGRAPH
CS IF(IGRAPH.GT.0)CALL INPLOT(NOW,F1,W)
IT=0
SI=0.D0
SI2=SI
SWGT=SI
SCHI=SI
SCALLS=SI
ENTRY VEGAS2(FXN,BCC,NDIM,NCALL,ITMX,NPRN,IGRAPH)
ND=NDMX
NG=1
IF(MDS.EQ.0) GO TO 2
NG=(NCALL*0.5)**(1./NDIM)
MDS=1
IF((2*NG-NDMX).LT.0) GO TO 2
MDS=-1
NPG=NG/NDMX+1
ND=NG/NPG
NG=NPG*ND
2 K=NG**NDIM
NPG=NCALL/K
IF(NPG.LT.2)NPG=2
CALLS=NPG*K
DXG=ONE/NG
DV2G=DXG**(2*NDIM)/NPG/NPG/(NPG-ONE)
XND=ND
NDM=ND-1
DXG=DXG*XND
XJAC=ONE
DO 3 J=1,NDIM
DX(J)=XU(J)-XL(J)
3 XJAC=XJAC*DX(J)
IF(ND.EQ.NDO) GO TO 8
RC=NDO/XND
DO 7 J=1,NDIM
K=0
XN=0.D0
DR=XN
I=K
4 K=K+1
DR=DR+ONE
XO=XN
XN=XI(K,J)
5 IF(RC.GT.DR) GO TO 4
I=I+1
DR=DR-RC
XIN(I)=XN-(XN-XO)*DR
IF(I.LT.NDM) GO TO 5
DO 6 I=1,NDM
6 XI(I,J)=XIN(I)
7 XI(ND,J)=ONE
NDO=ND
ACC=BCC
IF(NPRN.NE.0.AND.NPRN.NE.10)WRITE(NV,200)NDIM,CALLS,IT,ITMX
1,ACC,MDS,ND
8 CONTINUE
IF(NPRN.EQ.10)WRITE(NV,290)NDIM,CALLS,ITMX,ACC,MDS,ND
ENTRY VEGAS3(FXN,BCC,NDIM,NCALL,ITMX,NPRN,IGRAPH)
9 IT=IT+1
TI=0.D0
TSI=TI
CS IF(IGRAPH.GT.0)CALL REPLOT(NOW,F1,W)
DO 10 J=1,NDIM
KG(J)=1
DO 10 I=1,ND
NXI(I,J)=0
D(I,J)=TI
10 DI(I,J)=TI
11 FB=0.D0
F2B=FB
K=0
12 K=K+1
DO 121 J=1,NDIM
121 QRAN(J)=RANDM(0)
WGT=XJAC
DO 15 J=1,NDIM
XN=(KG(J)-QRAN(J))*DXG+ONE
IA(J)=XN
IAJ=IA(J)
IAJ1=IAJ-1
IF(IAJ.GT.1) GO TO 13
XO=XI(IAJ,J)
RC=(XN-IAJ)*XO
GO TO 14
13 XO=XI(IAJ,J)-XI(IAJ1,J)
RC=XI(IAJ1,J)+(XN-IAJ)*XO
14 X(J)=XL(J)+RC*DX(J)
15 WGT=WGT*XO*XND
F=FXN(X)*WGT
F1=F/CALLS
W=WGT/CALLS
CS IF(IGRAPH.GT.0)CALL XPLOT(NOW,F1,W)
F2=F**2
FB=FB+F
F2B=F2B+F2
DO 16 J=1,NDIM
IAJ=IA(J)
NXI(IAJ,J)=NXI(IAJ,J)+1
DI(IAJ,J)=DI(IAJ,J)+F/CALLS
16 IF(MDS.GE.0) D(IAJ,J)=D(IAJ,J)+F2
IF(K.LT.NPG) GO TO 12
F2B=F2B*NPG
F2B=SQRT(F2B)
F2B=(F2B-FB)*(F2B+FB)
TI=TI+FB
TSI=TSI+F2B
IF(MDS.GE.0) GO TO 18
DO 17 J=1,NDIM
IAJ=IA(J)
17 D(IAJ,J)=D(IAJ,J)+F2B
18 K=NDIM
19 KG(K)=MOD(KG(K),NG)+1
IF(KG(K).NE.1) GO TO 11
K=K-1
IF(K.GT.0) GO TO 19
TI=TI/CALLS
TSI=TSI*DV2G
TI2=TI*TI
WGT=TI2/TSI
SI=SI+TI*WGT
SI2=SI2+TI2
SWGT=SWGT+WGT
SCHI=SCHI+TI2*WGT
SCALLS=SCALLS+CALLS
AVGI=SI/SWGT
SD=SWGT*IT/SI2
CHI2A=0.D0
IF(IT.GT.1)CHI2A=SD*(SCHI/SWGT-AVGI*AVGI)/(IT-1)
SD=ONE/SD
SD=SQRT(SD)
IF(NPRN.EQ.0) GO TO 21
TSI=SQRT(TSI)
IF(NPRN.NE.10)WRITE(NV,201)IPR,IT,TI,TSI,AVGI,SD,CHI2A
RES=AVGI
IF(NPRN.EQ.10)WRITE(NV,203)IT,TI,TSI,AVGI,SD,CHI2A
IF(NPRN.GE.0) GO TO 21
DO 20 J=1,NDIM
WRITE(NV,202)J
20 WRITE(NV,204)(XI(I,J),DI(I,J),D(I,J),I=1,ND)
21 IF(ABS(SD/AVGI).LE.ABS(ACC).OR.IT.GE.ITMX)NOW=2
S1=AVGI
S2=SD
S3=TI
S4=TSI
CS IF(IGRAPH.GT.0)CALL PLOTIT(NOW,F1,W)
C DO 23 J=1,NDIM
C XO=D(1,J)
C XN=D(2,J)
C D(1,J)=(XO+XN)*0.5D0
C DT(J)=D(1,J)
C DO 22 I=2,NDM
C D(I,J)=XO+XN
C XO=XN
C XN=D(I+1,J)
C D(I,J)=(D(I,J)+XN)/3.D0
C22 DT(J)=DT(J)+D(I,J)
C D(ND,J)=(XN+XO)*0.5D0
C23 DT(J)=DT(J)+D(ND,J)
C-----THIS PART OF THE VEGAS-ALGORITHM IS UNSTABLE
C-----IT SHOULD BE REPLACED BY
DO 23 J=1,NDIM
DT(J)=0.D0
DO 23 I=1,ND
IF(NXI(I,J).GT.0)D(I,J)=D(I,J)/NXI(I,J)
23 DT(J)=DT(J)+D(I,J)
DO 28 J=1,NDIM
RC=0.D0
DO 24 I=1,ND
R(I)=0.D0
IF(D(I,J).LE.0.D0)GO TO 24
XO=DT(J)/D(I,J)
R(I)=((XO-ONE)/XO/LOG(XO))**ALPH
24 RC=RC+R(I)
RC=RC/XND
K=0
XN=0.D0
DR=XN
I=K
25 K=K+1
DR=DR+R(K)
XO=XN
XN=XI(K,J)
26 IF(RC.GT.DR) GO TO 25
I=I+1
DR=DR-RC
XIN(I)=XN-(XN-XO)*DR/R(K)
IF(I.LT.NDM) GO TO 26
DO 27 I=1,NDM
27 XI(I,J)=XIN(I)
28 XI(ND,J)=ONE
IF(IT.LT.ITMX.AND.ABS(ACC).LT.ABS(SD/AVGI))GO TO 9
200 FORMAT(35H0INPUT PARAMETERS FOR VEGAS NDIM=,I3
1,8H NCALL=,F8.0/28X,5H IT=,I5,8H ITMX =,I5/28X
2,6H ACC=,G9.3/28X,6H MDS=,I3,6H ND=,I4//)
290 FORMAT(13H0VEGAS NDIM=,I3,8H NCALL=,F8.0,8H ITMX =,I5
1,6H ACC=,G9.3,6H MDS=,I3,6H ND=,I4)
201 FORMAT(/I1,20HINTEGRATION BY VEGAS/13H0ITERATION NO,I3,
114H. INTEGRAL =,G14.8/20X,10HSTD DEV =,G10.4/
234H ACCUMULATED RESULTS. INTEGRAL =,G14.8/
324X,10HSTD DEV =,G10.4 / 24X,18HCHI**2 PER ITN =,G10.4)
202 FORMAT(14H0DATA FOR AXIS,I2 / 7X,1HX,7X,10H DELT I ,
12X,11H CONVCE ,11X,1HX,7X,10H DELT I ,2X,11H CONVCE
2,11X,1HX,7X,10H DELT I ,2X,11H CONVCE /)
204 FORMAT(1X,3G12.4,5X,3G12.4,5X,3G12.4)
203 FORMAT(1H ,I3,G20.8,G12.4,G20.8,G12.4,G12.4)
S1=AVGI
S2=SD
S3=CHI2A
RETURN
END
C----------------------------------------------------------------------
C A UNIVERSAL RANDOM NUMBER GENERATOR
DOUBLE PRECISION FUNCTION RANDM(IDMY)
C--RANDOM NUMBER GENERATOR
IMPLICIT REAL*8(A-H,O-Z)
REAL*4 UNIV
RANDM=DBLE(UNIV(1))
RETURN
END
C ---------------------------------------------------------------------
FUNCTION UNIV(IDUM)
C--FUNCTION FOR RANDOM NUMBER GENERATOR
REAL U(97)
COMMON /SET1/ U,C,CD,CM,I,J
UNIV=U(I)-U(J)
IF(UNIV.LT.0.) UNIV=UNIV+1.
U(I)=UNIV
I=I-1
IF(I.EQ.0) I=97
J=J-1
IF(J.EQ.0) J=97
C=C-CD
IF(C.LT.0.) C=C+CM
UNIV=UNIV-C
IF(UNIV.LT.0.) UNIV=UNIV+1
RETURN
END
C----------------------------------------------------------------------
C INITIALIZING THE RANDOM NUMBER GENERATOR
C TO INITIALIZE CALL RSTART(12,34,56,78)
SUBROUTINE RSTART(I,J,K,L)
C--INITIALIZATION ROUTINE FOR RANDOM NUMBER GENERATOR
REAL U(97)
COMMON /SET1/ U,C,CD,CM,ISTART,JSTART
IF ((I.LT.0).OR.(I.GT.178 )) STOP 'FIRST SEED .LT.0 OR .GT.178'
IF ((J.LT.0).OR.(J.GT.178 )) STOP 'SECOND SEED .LT.0 OR .GT.178'
IF ((K.LT.0).OR.(K.GT.178 )) STOP 'THIRD SEED .LT.0 OR .GT.178'
IF ((L.LT.0).OR.(L.GT.168 )) STOP 'FOURTH SEED .LT.0 OR .GT.168'
IF ( (I.EQ.1).AND.(J.EQ.1).AND.(K.EQ.1) ) STOP
& 'FIRST, SECOND AND THIRD SEEDS ARE ALL EQUAL TO 1'
ISTART=97
JSTART=33
IDUM=I
JDUM=J
KDUM=K
LDUM=L
DO 2 II=1,97
S=0.
T=.5
DO 3 JJ=1,24
M=MOD(MOD(IDUM*JDUM,179)*K,179)
IDUM=JDUM
JDUM=KDUM
KDUM=M
LDUM=MOD(53*LDUM+1,169)
IF(MOD(LDUM*M,64).GE.32) S=S+T
3 T=.5*T
2 U(II)=S
C=362436./16777216.
CD=7654321./16777216.
CM=16777213./16777216.
RETURN
END
C%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
C *****************************************************************
C ************* SUBROUTINE FOR THE SUSY COUPLINGS *****************
C *****************************************************************
SUBROUTINE SUSYCP(TGBET)
IMPLICIT DOUBLE PRECISION (A-H,O-Z)
DOUBLE PRECISION LA1,LA2,LA3,LA4,LA5,LA6,LA7,LA3T
COMMON/MASSES/AMS,AMC,AMB,AMT
COMMON/HMASS/AMSM,AMA,AML,AMH,AMCH
COMMON/HSELF/LA1,LA2,LA3,LA4,LA5,LA6,LA7
COMMON/BREAK/AMSQ,AMUR,AU,AD,AMU,AM2,ASL
COMMON/PARAM/S,AM,AMQ,SCALE,EPS,GF,GEVPB
COMMON/VMASSES/AMW,AMZ
COMMON/COUP/GAT,GAB,GLT,GLB,GHT,GHB,GZAH,GZAL,
. GHHH,GLLL,GHLL,GLHH,GHAA,GLAA,GLVV,GHVV,
. GLPM,GHPM,B,A
PI=4*DATAN(1D0)
V=1.D0/DSQRT(DSQRT(2.D0)*GF)
BET=DATAN(TGBET)
SB = DSIN(BET)
CB = DCOS(BET)
CALL SUBH(AMA,TGBET,AMSQ,AMUR,AMT,AU,AD,AMU,AML,AMH,AMCH,SA,CA,
. TBA)
LA3T=LA3+LA4+LA5
AMA2=AMA**2
AML2=AML**2
AMH2=AMH**2
AMP2=AMCH**2
SBMA = SB*CA-CB*SA
CBMA = CB*CA+SB*SA
SBPA = SB*CA+CB*SA
CBPA = CB*CA-SB*SA
S2A = 2*SA*CA
C2A = CA**2-SA**2
S2B = 2*SB*CB
C2B = CB**2-SB**2
GLZZ = 1/V/2*AML2*SBMA
GHZZ = 1/V/2*AMH2*CBMA
GLWW = 2*GLZZ
GHWW = 2*GHZZ
GLAZ = 1/V*(AML2-AMA2)*CBMA
GHAZ = -1/V*(AMH2-AMA2)*SBMA
GLPW = -1/V*(AMP2-AML2)*CBMA
GLMW = GLPW
GHPW = 1/V*(AMP2-AMH2)*SBMA
GHMW = GHPW
GAPW = 1/V*(AMP2-AMA2)
GAMW = -GAPW
GHHH = V/2*(LA1*CA**3*CB + LA2*SA**3*SB + LA3T*SA*CA*SBPA
. + LA6*CA**2*(3*SA*CB+CA*SB) + LA7*SA**2*(3*CA*SB+SA*CB))
GLLL = -V/2*(LA1*SA**3*CB - LA2*CA**3*SB + LA3T*SA*CA*CBPA
. - LA6*SA**2*(3*CA*CB-SA*SB) + LA7*CA**2*(3*SA*SB-CA*CB))
GLHH = -3*V/2*(LA1*CA**2*CB*SA - LA2*SA**2*SB*CA
. + LA3T*(SA**3*CB-CA**3*SB+2*SBMA/3)
. - LA6*CA*(CB*C2A-SA*SBPA) - LA7*SA*(C2A*SB+CA*SBPA))
GHLL = 3*V/2*(LA1*SA**2*CB*CA + LA2*CA**2*SB*SA
. + LA3T*(SA**3*SB+CA**3*CB-2*CBMA/3)
. - LA6*SA*(CB*C2A+CA*CBPA) + LA7*CA*(C2A*SB+SA*CBPA))
GLAA = -V/2*(LA1*SB**2*CB*SA - LA2*CB**2*SB*CA
. - LA3T*(SB**3*CA-CB**3*SA) + 2*LA5*SBMA
. - LA6*SB*(CB*SBPA+SA*C2B) - LA7*CB*(C2B*CA-SB*SBPA))
GHAA = V/2*(LA1*SB**2*CB*CA + LA2*CB**2*SB*SA
. + LA3T*(SB**3*SA+CB**3*CA) - 2*LA5*CBMA
. - LA6*SB*(CB*CBPA+CA*C2B) + LA7*CB*(SB*CBPA+SA*C2B))
GLPM = 2*GLAA + V*(LA5 - LA4)*SBMA
GHPM = 2*GHAA + V*(LA5 - LA4)*CBMA
GLZZ = 2*GLZZ
GHZZ = 2*GHZZ
GLLL = 6*GLLL
GHHH = 6*GHHH
GLHH = 2*GLHH
GHLL = 2*GHLL
GLAA = 2*GLAA
GHAA = 2*GHAA
XNORM = AMZ**2/V
GLLL = GLLL/XNORM
GHLL = GHLL/XNORM
GLHH = GLHH/XNORM
GHHH = GHHH/XNORM
GHAA = GHAA/XNORM
GLAA = GLAA/XNORM
GLPM = GLPM/XNORM
GHPM = GHPM/XNORM
GAT=1.D0/TGBET
GAB=TGBET
GLT=CA/SB
GLB=-SA/CB
GHT=SA/SB
GHB=CA/CB
GZAL=-CBMA
GZAH=SBMA
GLVV=SBMA
GHVV=CBMA
B=BET
A=DATAN(SA/CA)
IF(CA.LT.0D0)THEN
IF(SA.LT.0D0)THEN
A = A-PI
ELSE
A = A+PI
ENDIF
ENDIF
RETURN
END
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C THIS PROGRAM COMPUTES THE RENORMALIZATION GROUP IMPROVED
C VALUES OF HIGGS MASSES AND COUPLINGS IN THE MSSM.
C
C INPUT: MA,TANB = TAN(BETA),MQ,MUR,MTOP,AU,AD,MU.
C
C ALL MASSES IN GEV UNITS. MA IS THE CP-ODD HIGGS MASS,
C MTOP IS THE PHYSICAL TOP MASS, MQ AND MUR ARE THE SOFT
C SUPERSYMMETRY BREAKING MASS PARAMETERS OF LEFT HANDED
C AND RIGHT HANDED STOPS RESPECTIVELY, AU AND AD ARE THE
C STOP AND SBOTTOM TRILINEAR SOFT BREAKING TERMS,
C RESPECTIVELY, AND MU IS THE SUPERSYMMETRIC
C HIGGS MASS PARAMETER. WE USE THE CONVENTIONS FROM
C THE PHYSICS REPORT OF HABER AND KANE: LEFT RIGHT
C STOP MIXING TERM PROPORTIONAL TO (AU - MU/TANB).
C
C WE USE AS INPUT TANB DEFINED AT THE SCALE MTOP.
C
C OUTPUT: MH,HM,MHCH, SA = SIN(ALPHA), CA= COS(ALPHA), TANBA
C
C WHERE MH AND HM ARE THE LIGHTEST AND HEAVIEST CP-EVEN
C HIGGS MASSES, MHCH IS THE CHARGED HIGGS MASS AND
C ALPHA IS THE HIGGS MIXING ANGLE.
C
C TANBA IS THE ANGLE TANB AT THE CP-ODD HIGGS MASS SCALE.
C
C RANGE OF VALIDITY:
C
C (STOP1**2 - STOP2**2)/(STOP2**2 + STOP1**2) < 0.5
C (SBOT1**2 - SBOT2**2)/(SBOT2**2 + SBOT2**2) < 0.5
C
C WHERE STOP1, STOP2, SBOT1 AND SBOT2 ARE THE STOP AND
C ARE THE SBOTTOM MASS EIGENVALUES, RESPECTIVELY. THIS
C RANGE AUTOMATICALLY EXCLUDES THE EXISTENCE OF TACHYONS.
C
C
C FOR THE CHARGED HIGGS MASS COMPUTATION, THE METHOD IS
C VALID IF
C
C 2 * |MB * AD* TANB| < M_SUSY**2, 2 * |MTOP * AU| < M_SUSY**2
C
C 2 * |MB * MU * TANB| < M_SUSY**2, 2 * |MTOP * MU| < M_SUSY**2
C
C WHERE M_SUSY**2 IS THE AVERAGE OF THE SQUARED STOP MASS
C EIGENVALUES, M_SUSY**2 = (STOP1**2 + STOP2**2)/2. THE SBOTTOM
C MASSES HAVE BEEN ASSUMED TO BE OF ORDER OF THE STOP ONES.
C
C M_SUSY**2 = (MQ**2 + MUR**2)*0.5 + MTOP**2
C
C PROGRAM BASED ON THE WORK BY M. CARENA, J.R. ESPINOSA,
C M. QUIROS AND C.E.M. WAGNER, CERN-PREPRINT CERN-TH/95-45.
C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
SUBROUTINE SUBH (MA,TANB,MQ,MUR,MTOP,AU,AD,MU,MH,HM,
* MHCH,SA,CA,TANBA)
IMPLICIT REAL*8(A-H,L,M,O-Z)
COMMON/PARAM/S,AM,AMQ,SCALE,EPS,GF,GEVPB
COMMON/VMASSES/AMW,AMZ
COMMON/HSELF/LAMBDA1,LAMBDA2,LAMBDA3,LAMBDA4,LAMBDA5,
. LAMBDA6,LAMBDA7
C MZ = 91.18
C ALPHA1 = 0.0101
C ALPHA2 = 0.0337
C ALPHA3Z = 0.12
C V = 174.1
C PI = 3.14159
TANBA = TANB
TANBT = TANB
C MBOTTOM(MTOP) = 3. GEV
C MB = 3.
C ALPHA3 = ALPHA3Z/(1. +(11. - 10./3.)/4./PI*ALPHA3Z*
C *LOG(MTOP**2/MZ**2))
C RMTOP= RUNNING TOP QUARK MASS
C RMTOP = MTOP/(1.+4.*ALPHA3/3./PI)
C>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
MB = RUNM(MTOP,5)
PI = 4*DATAN(1D0)
MZ = AMZ
V = 1/DSQRT(2*DSQRT(2D0)*GF)
CW = AMW**2/AMZ**2
SW = 1-CW
ALPHA2 = (2*AMW/V/DSQRT(2D0))**2/4/PI
ALPHA1 = ALPHA2*SW/CW
ALPHA3Z = ALPHAS(AMZ,2)
ALPHA3 = ALPHAS(MTOP,2)
RMTOP = RUNM(MTOP,6)
C>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
C RMTOP=MTOP
C>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
MS = ((MQ**2 + MUR**2)/2. + MTOP**2)**.5
T = LOG(MS**2/MTOP**2)
SINB = TANB/((1. + TANB**2)**.5)
COSB = SINB/TANB
C IF(MA.LE.MTOP) TANBA = TANBT
IF(MA.GT.MTOP)
*TANBA = TANBT*(1.-3./32./PI**2*
*(RMTOP**2/V**2/SINB**2-MB**2/V**2/COSB**2)*
*LOG(MA**2/MTOP**2))
SINBT = TANBT/((1. + TANBT**2)**.5)
COSBT = 1./((1. + TANBT**2)**.5)
COS2BT = (TANBT**2 - 1.)/(TANBT**2 + 1.)
G1 = (ALPHA1*4.*PI)**.5
G2 = (ALPHA2*4.*PI)**.5
G3 = (ALPHA3*4.*PI)**.5
HU = RMTOP/V/SINBT
HD = MB/V/COSBT
C>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
C G3=0
C HU=0
C HD=0
C>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
XAU = (2.*AU**2/MS**2)*(1. - AU**2/12./MS**2)
XAD = (2.*AD**2/MS**2)*(1. - AD**2/12./MS**2)
AUD = (-6.*MU**2/MS**2 - ( MU**2- AD*AU)**2/MS**4.
*+ 3.*(AU + AD)**2/MS**2)/6.
LAMBDA1 = ((G1**2 + G2**2)/4.)*(1.-3.*HD**2*T/8./PI**2)
*+(3.*HD**4/8./PI**2) * (T + XAD/2. + (3.*HD**2/2. + HU**2/2.
*- 8.*G3**2) * (XAD*T + T**2)/16./PI**2)
*-(3.*HU**4* MU**4/96./PI**2/MS**4) * (1+ (9.*HU**2 -5.* HD**2
*- 16.*G3**2) *T/16./PI**2)
LAMBDA2 = ((G1**2 + G2**2)/4.)*(1.-3.*HU**2*T/8./PI**2)
*+(3.*HU**4/8./PI**2) * (T + XAU/2. + (3.*HU**2/2. + HD**2/2.
*- 8.*G3**2) * (XAU*T + T**2)/16./PI**2)
*-(3.*HD**4* MU**4/96./PI**2/MS**4) * (1+ (9.*HD**2 -5.* HU**2
*- 16.*G3**2) *T/16./PI**2)
LAMBDA3 = ((G2**2 - G1**2)/4.)*(1.-3.*
*(HU**2 + HD**2)*T/16./PI**2)
*+(6.*HU**2*HD**2/16./PI**2) * (T + AUD/2. + (HU**2 + HD**2
*- 8.*G3**2) * (AUD*T + T**2)/16./PI**2)
*+(3.*HU**4/96./PI**2) * (3.*MU**2/MS**2 - MU**2*AU**2/
*MS**4)* (1.+ (6.*HU**2 -2.* HD**2/2.
*- 16.*G3**2) *T/16./PI**2)
*+(3.*HD**4/96./PI**2) * (3.*MU**2/MS**2 - MU**2*AD**2/
*MS**4)*(1.+ (6.*HD**2 -2.* HU**2
*- 16.*G3**2) *T/16./PI**2)
LAMBDA4 = (- G2**2/2.)*(1.-3.*(HU**2 + HD**2)*T/16./PI**2)
*-(6.*HU**2*HD**2/16./PI**2) * (T + AUD/2. + (HU**2 + HD**2
*- 8.*G3**2) * (AUD*T + T**2)/16./PI**2)
*+(3.*HU**4/96./PI**2) * (3.*MU**2/MS**2 - MU**2*AU**2/
*MS**4)*
*(1+ (6.*HU**2 -2.* HD**2
*- 16.*G3**2) *T/16./PI**2)
*+(3.*HD**4/96./PI**2) * (3.*MU**2/MS**2 - MU**2*AD**2/
*MS**4)*
*(1+ (6.*HD**2 -2.* HU**2/2.
*- 16.*G3**2) *T/16./PI**2)
LAMBDA5 = -(3.*HU**4* MU**2*AU**2/96./PI**2/MS**4) *
* (1- (2.*HD**2 -6.* HU**2 + 16.*G3**2) *T/16./PI**2)
*-(3.*HD**4* MU**2*AD**2/96./PI**2/MS**4) *
* (1- (2.*HU**2 -6.* HD**2 + 16.*G3**2) *T/16./PI**2)
LAMBDA6 = (3.*HU**4* MU**3*AU/96./PI**2/MS**4) *
* (1- (7.*HD**2/2. -15.* HU**2/2. + 16.*G3**2) *T/16./PI**2)
*+(3.*HD**4* MU *(AD**3/MS**3 - 6.*AD/MS )/96./PI**2/MS) *
* (1- (HU**2/2. -9.* HD**2/2. + 16.*G3**2) *T/16./PI**2)
LAMBDA7 = (3.*HD**4* MU**3*AD/96./PI**2/MS**4) *
* (1- (7.*HU**2/2. -15.* HD**2/2. + 16.*G3**2) *T/16./PI**2)
*+(3.*HU**4* MU *(AU**3/MS**3 - 6.*AU/MS )/96./PI**2/MS) *
* (1- (HD**2/2. -9.* HU**2/2. + 16.*G3**2) *T/16./PI**2)
TRM2 = MA**2 + 2.*V**2* (LAMBDA1* COSBT**2 +
*2.* LAMBDA6*SINBT*COSBT
*+ LAMBDA5*SINBT**2 + LAMBDA2* SINBT**2 + 2.* LAMBDA7*SINBT*COSBT
*+ LAMBDA5*COSBT**2)
DETM2 = 4.*V**4*(-(SINBT*COSBT*(LAMBDA3 + LAMBDA4) +
*LAMBDA6*COSBT**2
*+ LAMBDA7* SINBT**2)**2 + (LAMBDA1* COSBT**2 +
*2.* LAMBDA6* COSBT*SINBT
*+ LAMBDA5*SINBT**2)*(LAMBDA2* SINBT**2 +2.* LAMBDA7* COSBT*SINBT
*+ LAMBDA5*COSBT**2)) + MA**2*2.*V**2 *
*((LAMBDA1* COSBT**2 +2.*
*LAMBDA6* COSBT*SINBT + LAMBDA5*SINBT**2)*COSBT**2 +
*(LAMBDA2* SINBT**2 +2.* LAMBDA7* COSBT*SINBT + LAMBDA5*COSBT**2)
**SINBT**2
* +2.*SINBT*COSBT* (SINBT*COSBT*(LAMBDA3
* + LAMBDA4) + LAMBDA6*COSBT**2
*+ LAMBDA7* SINBT**2))
MH2 = (TRM2 - (TRM2**2 - 4.* DETM2)**.5)/2.
HM2 = (TRM2 + (TRM2**2 - 4.* DETM2)**.5)/2.
HM = HM2**.5
MH = MH2**.5
MHCH2 = MA**2 + (LAMBDA5 - LAMBDA4)* V**2
MHCH = MHCH2**.5
MHCH = MHCH2**.5
SINALPHA = SQRT(((TRM2**2 - 4.* DETM2)**.5) -
* ((2.*V**2*(LAMBDA1* COSBT**2 + 2.*
*LAMBDA6* COSBT*SINBT
*+ LAMBDA5*SINBT**2) + MA**2*SINBT**2)
*- (2.*V**2*(LAMBDA2* SINBT**2 +2.* LAMBDA7* COSBT*SINBT
*+ LAMBDA5*COSBT**2) + MA**2*COSBT**2)))/
*SQRT(((TRM2**2 - 4.* DETM2)**.5))/2.**.5
COSALPHA = (2.*(2.*V**2*(SINBT*COSBT*(LAMBDA3 + LAMBDA4) +
*LAMBDA6*COSBT**2 + LAMBDA7* SINBT**2) -
*MA**2*SINBT*COSBT))/2.**.5/
*SQRT(((TRM2**2 - 4.* DETM2)**.5)*
*(((TRM2**2 - 4.* DETM2)**.5) -
* ((2.*V**2*(LAMBDA1* COSBT**2 + 2.*
*LAMBDA6* COSBT*SINBT
*+ LAMBDA5*SINBT**2) + MA**2*SINBT**2)
*- (2.*V**2*(LAMBDA2* SINBT**2 +2.* LAMBDA7* COSBT*SINBT
*+ LAMBDA5*COSBT**2) + MA**2*COSBT**2))))
SA = -SINALPHA
CA = -COSALPHA
2242 RETURN
END
SUBROUTINE AMHAMA (ICASE,MH,MA,TANB,MQ,MUR,MTOP,AU,AD,MU)
C--CALCULATION OF PSEUDOSCALAR HIGGS MASS FROM HIGGS MASS MH
C--ICASE=0: MH=PSEUDOSCALAR MASS
C--ICASE=1: MH=LIGHT SCALAR MASS
C--ICASE=2: MH=HEAVY SCALAR MASS
C--ICASE=3: MH=CHARGED HIGGS MASS
IMPLICIT REAL*8(A-H,L,M,O-Z)
COMMON/PARAM/S,AM,AMQ,SCALE,EPS,GF,GEVPB
COMMON/VMASSES/AMW,AMZ
C MZ = 91.18
C ALPHA1 = 0.0101
C ALPHA2 = 0.0337
C ALPHA3Z = 0.12
C V = 174.1
C PI = 3.14159
TANBI = TANB
TANBA = TANB
TANBT = TANB
C MBOTTOM(MTOP) = 3. GEV
C MB = 3.
C ALPHA3 = ALPHA3Z/(1. +(11. - 10./3.)/4./PI*ALPHA3Z*
C *LOG(MTOP**2./MZ**2.))
C RMTOP= RUNNING TOP QUARK MASS
C RMTOP = MTOP/(1.+4.*ALPHA3/3./PI)
C>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
MB = RUNM(MTOP,5)
PI = 4*DATAN(1D0)
MZ = AMZ
V = 1/DSQRT(2*DSQRT(2D0)*GF)
CW = AMW**2/AMZ**2
SW = 1-CW
ALPHA2 = (2*AMW/V/DSQRT(2D0))**2/4/PI
ALPHA1 = ALPHA2*SW/CW
ALPHA3Z = ALPHAS(AMZ,2)
ALPHA3 = ALPHAS(MTOP,2)
RMTOP = RUNM(MTOP,6)
C>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
C>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
C RMTOP=MTOP
C>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
MS = ((MQ**2. + MUR**2.)/2. + MTOP**2.)**.5
T = LOG(MS**2./MTOP**2.)
SINB = TANBI/((1. + TANBI**2.)**.5)
COSB = SINB/TANBI
C IF(MA.LE.MTOP) TANBT = TANBI
IF(MA.GT.MTOP)
*TANBA = TANBT*(1.-3./32./PI**2.*
*(RMTOP**2./V**2./SINB**2.-MB**2./V**2./COSB**2.)*
*LOG(MA**2./MTOP**2.))
SINBT = TANBT/((1. + TANBT**2.)**.5)
COSBT = 1./((1. + TANBT**2.)**.5)
COS2BT = (TANBT**2. - 1.)/(TANBT**2. + 1.)
G1 = (ALPHA1*4.*PI)**.5
G2 = (ALPHA2*4.*PI)**.5
G3 = (ALPHA3*4.*PI)**.5
HU = RMTOP/V/SINBT
HD = MB/V/COSBT
C>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
C G3=0
C HU=0
C HD=0
C>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
XAU = (2.*AU**2./MS**2.)*(1. - AU**2./12./MS**2.)
XAD = (2.*AD**2./MS**2.)*(1. - AD**2./12./MS**2.)
AUD = (-6.*MU**2/MS**2. - ( MU**2- AD*AU)**2./MS**4.
*+ 3.*(AU + AD)**2./MS**2.)/6.
LAMBDA1 = ((G1**2. + G2**2.)/4.)*(1.-3.*HD**2.*T/8./PI**2.)
*+(3.*HD**4./8./PI**2.) * (T + XAD/2. + (3.*HD**2./2. + HU**2./2.
*- 8.*G3**2.) * (XAD*T + T**2.)/16./PI**2.)
*-(3.*HU**4.* MU**4/96./PI**2./MS**4.) * (1+ (9.*HU**2. -5.* HD**2.
*- 16.*G3**2.) *T/16./PI**2.)
LAMBDA2 = ((G1**2. + G2**2.)/4.)*(1.-3.*HU**2.*T/8./PI**2.)
*+(3.*HU**4./8./PI**2.) * (T + XAU/2. + (3.*HU**2./2. + HD**2./2.
*- 8.*G3**2.) * (XAU*T + T**2.)/16./PI**2.)
*-(3.*HD**4.* MU**4/96./PI**2./MS**4.) * (1+ (9.*HD**2. -5.* HU**2.
*- 16.*G3**2.) *T/16./PI**2.)
LAMBDA3 = ((G2**2. - G1**2.)/4.)*(1.-3.*
*(HU**2. + HD**2.)*T/16./PI**2.)
*+(6.*HU**2.*HD**2./16./PI**2.) * (T + AUD/2. + (HU**2. + HD**2.
*- 8.*G3**2.) * (AUD*T + T**2.)/16./PI**2.)
*+(3.*HU**4./96./PI**2.) * (3.*MU**2/MS**2. - MU**2*AU**2./
*MS**4.)* (1.+ (6.*HU**2. -2.* HD**2./2.
*- 16.*G3**2.) *T/16./PI**2.)
*+(3.*HD**4./96./PI**2.) * (3.*MU**2/MS**2. - MU**2*AD**2./
*MS**4.)*(1.+ (6.*HD**2. -2.* HU**2./2.
*- 16.*G3**2.) *T/16./PI**2.)
LAMBDA4 = (- G2**2./2.)*(1.-3.*(HU**2. + HD**2.)*T/16./PI**2.)
*-(6.*HU**2.*HD**2./16./PI**2.) * (T + AUD/2. + (HU**2. + HD**2.
*- 8.*G3**2.) * (AUD*T + T**2.)/16./PI**2.)
*+(3.*HU**4./96./PI**2.) * (3.*MU**2/MS**2. - MU**2*AU**2./
*MS**4.)*
*(1+ (6.*HU**2. -2.* HD**2./2.
*- 16.*G3**2.) *T/16./PI**2.)
*+(3.*HD**4./96./PI**2.) * (3.*MU**2/MS**2. - MU**2*AD**2./
*MS**4.)*
*(1+ (6.*HD**2. -2.* HU**2./2.
*- 16.*G3**2.) *T/16./PI**2.)
LAMBDA5 = -(3.*HU**4.* MU**2*AU**2./96./PI**2./MS**4.) *
* (1- (2.*HD**2. -6.* HU**2. + 16.*G3**2.) *T/16./PI**2.)
*-(3.*HD**4.* MU**2*AD**2./96./PI**2./MS**4.) *
* (1- (2.*HU**2. -6.* HD**2. + 16.*G3**2.) *T/16./PI**2.)
LAMBDA6 = (3.*HU**4.* MU**3*AU/96./PI**2./MS**4.) *
* (1- (7.*HD**2./2. -15.* HU**2./2. + 16.*G3**2.) *T/16./PI**2.)
*+(3.*HD**4.* MU *(AD**3./MS**3. - 6.*AD/MS )/96./PI**2./MS) *
* (1- (HU**2./2. -9.* HD**2./2. + 16.*G3**2.) *T/16./PI**2.)
LAMBDA7 = (3.*HD**4.* MU**3*AD/96./PI**2./MS**4.) *
* (1- (7.*HU**2./2. -15.* HD**2./2. + 16.*G3**2.) *T/16./PI**2.)
*+(3.*HU**4.* MU *(AU**3./MS**3. - 6.*AU/MS )/96./PI**2./MS) *
* (1- (HD**2./2. -9.* HU**2./2. + 16.*G3**2.) *T/16./PI**2.)
DEL1 = 2.*V**2.* (LAMBDA1* COSBT**2 +
*2.* LAMBDA6*SINBT*COSBT
*+ LAMBDA5*SINBT**2. + LAMBDA2* SINBT**2 + 2.* LAMBDA7*SINBT*COSBT
*+ LAMBDA5*COSBT**2.)
DEL2 = 4.*V**4.*(-(SINBT*COSBT*(LAMBDA3 + LAMBDA4) +
*LAMBDA6*COSBT**2.
*+ LAMBDA7* SINBT**2.)**2. + (LAMBDA1* COSBT**2. +
*2.* LAMBDA6* COSBT*SINBT
*+ LAMBDA5*SINBT**2.)*(LAMBDA2* SINBT**2. +2.* LAMBDA7* COSBT*SINBT
*+ LAMBDA5*COSBT**2.))
DELA = 2.*V**2. *
*((LAMBDA1* COSBT**2. +2.*
*LAMBDA6* COSBT*SINBT + LAMBDA5*SINBT**2.)*COSBT**2. +
*(LAMBDA2* SINBT**2. +2.* LAMBDA7* COSBT*SINBT + LAMBDA5*COSBT**2.)
**SINBT**2.
* +2.*SINBT*COSBT* (SINBT*COSBT*(LAMBDA3
* + LAMBDA4) + LAMBDA6*COSBT**2.
*+ LAMBDA7* SINBT**2.))
IF(ICASE.EQ.1.OR.ICASE.EQ.2)THEN
XX = (MH**4 - MH**2*DEL1 + DEL2)/(MH**2 - DELA)
IF(XX.GT.0D0)THEN
MA = DSQRT(XX)
ELSE
MA = -1.D8
ENDIF
ELSEIF(ICASE.EQ.3)THEN
XX = MH**2 - (LAMBDA5 - LAMBDA4)* V**2.
IF(XX.GT.0D0)THEN
MA = DSQRT(MH**2 - (LAMBDA5 - LAMBDA4)* V**2.)
ELSE
MA = -1.D8
ENDIF
ELSE
MA = MH
ENDIF
RETURN
END
DOUBLE PRECISION FUNCTION RUNM(Q,NF)
PARAMETER (NN=6)
IMPLICIT DOUBLE PRECISION (A-H,O-Z)
PARAMETER (ZETA3 = 1.202056903159594D0)
DIMENSION AM(NN),YMSB(NN)
COMMON/ALS/XLAMBDA,AMCA,AMBA,AMTA,N0A
COMMON/MASSES/AMS,AMC,AMB,AMT
COMMON/STRANGE/AMSB
COMMON/RUN/XMSB,XMHAT,XKFAC
SAVE ISTRANGE
B0(NF)=(33.D0-2.D0*NF)/12D0
B1(NF) = (102D0-38D0/3D0*NF)/16D0
B2(NF) = (2857D0/2D0-5033D0/18D0*NF+325D0/54D0*NF**2)/64D0
G0(NF) = 1D0
G1(NF) = (202D0/3D0-20D0/9D0*NF)/16D0
G2(NF) = (1249D0-(2216D0/27D0+160D0/3D0*ZETA3)*NF
. - 140D0/81D0*NF**2)/64D0
C1(NF) = G1(NF)/B0(NF) - B1(NF)*G0(NF)/B0(NF)**2
C2(NF) = ((G1(NF)/B0(NF) - B1(NF)*G0(NF)/B0(NF)**2)**2
. + G2(NF)/B0(NF) + B1(NF)**2*G0(NF)/B0(NF)**3
. - B1(NF)*G1(NF)/B0(NF)**2 - B2(NF)*G0(NF)/B0(NF)**2)/2D0
TRAN(X,XK)=1D0+4D0/3D0*ALPHAS(X,2)/PI+XK*(ALPHAS(X,2)/PI)**2
CQ(X,NF)=(2D0*B0(NF)*X)**(G0(NF)/B0(NF))
. *(1D0+C1(NF)*X+C2(NF)*X**2)
DATA ISTRANGE/0/
PI=4D0*DATAN(1D0)
ACC = 1.D-8
AM(1) = 0
AM(2) = 0
C--------------------------------------------
IMSBAR = 0
NNLO = 0
IF(IMSBAR.EQ.1)THEN
IF(ISTRANGE.EQ.0)THEN
C--STRANGE POLE MASS FROM MSBAR-MASS AT 1 GEV
AMSD = XLAMBDA
AMSU = 1.D8
123 AMS = (AMSU+AMSD)/2
AM(3) = AMS
XMSB = AMS/CQ(ALPHAS(AMS,2)/PI,3)
. *CQ(ALPHAS(1.D0,2)/PI,3)/TRAN(AMS,0D0)
DD = (XMSB-AMSB)/AMSB
IF(DABS(DD).GE.ACC)THEN
IF(DD.LE.0.D0)THEN
AMSD = AM(3)
ELSE
AMSU = AM(3)
ENDIF
GOTO 123
ENDIF
ISTRANGE=1
ENDIF
AM(3) = AMSB
ELSE
AMS=AMSB
AM(3) = AMS
ENDIF
C--------------------------------------------
AM(3) = AMSB
AM(4) = AMC
AM(5) = AMB
AM(6) = AMT
XK = 16.11D0
DO 1 I=1,NF-1
XK = XK - 1.04D0*(1.D0-AM(I)/AM(NF))
1 CONTINUE
IF(NF.GE.4)THEN
XMSB = AM(NF)/TRAN(AM(NF),0D0)
XMHAT = XMSB/CQ(ALPHAS(AM(NF),2)/PI,NF)
ELSE
XMSB = 0
XMHAT = 0
ENDIF
YMSB(3) = AMSB
IF(NF.EQ.3)THEN
YMSB(4) = YMSB(3)*CQ(ALPHAS(AM(4),2)/PI,3)/
. CQ(ALPHAS(1.D0,2)/PI,3)
YMSB(5) = YMSB(4)*CQ(ALPHAS(AM(5),2)/PI,4)/
. CQ(ALPHAS(AM(4),2)/PI,4)
YMSB(6) = YMSB(5)*CQ(ALPHAS(AM(6),2)/PI,5)/
. CQ(ALPHAS(AM(5),2)/PI,5)
ELSEIF(NF.EQ.4)THEN
YMSB(4) = XMSB
YMSB(5) = YMSB(4)*CQ(ALPHAS(AM(5),2)/PI,4)/
. CQ(ALPHAS(AM(4),2)/PI,4)
YMSB(6) = YMSB(5)*CQ(ALPHAS(AM(6),2)/PI,5)/
. CQ(ALPHAS(AM(5),2)/PI,5)
ELSEIF(NF.EQ.5)THEN
YMSB(5) = XMSB
YMSB(4) = YMSB(5)*CQ(ALPHAS(AM(4),2)/PI,4)/
. CQ(ALPHAS(AM(5),2)/PI,4)
YMSB(6) = YMSB(5)*CQ(ALPHAS(AM(6),2)/PI,5)/
. CQ(ALPHAS(AM(5),2)/PI,5)
ELSEIF(NF.EQ.6)THEN
YMSB(6) = XMSB
YMSB(5) = YMSB(6)*CQ(ALPHAS(AM(5),2)/PI,5)/
. CQ(ALPHAS(AM(6),2)/PI,5)
YMSB(4) = YMSB(5)*CQ(ALPHAS(AM(4),2)/PI,4)/
. CQ(ALPHAS(AM(5),2)/PI,4)
ENDIF
IF(Q.LT.AMC)THEN
N0=3
Q0 = 1.D0
ELSEIF(Q.LE.AMB)THEN
N0=4
Q0 = AMC
ELSEIF(Q.LE.AMT)THEN
N0=5
Q0 = AMB
ELSE
N0=6
Q0 = AMT
ENDIF
IF(NNLO.EQ.1.AND.NF.GT.3)THEN
XKFAC = TRAN(AM(NF),0D0)/TRAN(AM(NF),XK)
ELSE
XKFAC = 1D0
ENDIF
RUNM = YMSB(N0)*CQ(ALPHAS(Q,2)/PI,N0)/
. CQ(ALPHAS(Q0,2)/PI,N0)
. * XKFAC
RETURN
END
C%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
C SUBROUTINE STRUC(X,Q,PDF)
C--PARTON DENSITIES
C IMPLICIT DOUBLE PRECISION (A-H,O-Z)
C DIMENSION PDF(-6:6), VALUE(20)
C CHARACTER*20 PARM(20)
C COMMON/PDFLIB/NGROUP,NSET
C QQ=Q
C PARM(1)='NPTYPE'
C PARM(2)='NGROUP'
C PARM(3)='NSET'
C VALUE(1)=1.D0
C VALUE(2)=DFLOAT(NGROUP)
C VALUE(3)=DFLOAT(NSET)
C CALL PDFSET(PARM,VALUE)
C CALL PFTOPDG(X,QQ,PDF)
C RETURN
C END
subroutine struc(x,q,pdf)
implicit double precision (a-h,o-z)
dimension pdf(-6:6), value(20)
character*20 parm(20)
common/pdflib/ngroup,nset
if(ngroup.gt.0)then
parm(1)='nptype'
parm(2)='ngroup'
parm(3)='nset'
value(1)=1.d0
value(2)=dfloat(ngroup)
value(3)=dfloat(nset)
call pdfset(parm,value)
call pftopdg(x,q,pdf)
elseif(ngroup.eq.-1)then
call SetCtq6(nset)
pdf(6) = 0
pdf(-6) = 0
do i=-5,5
j = i
if(i.eq.1)j=2
if(i.eq.2)j=1
if(i.eq.-1)j=-2
if(i.eq.-2)j=-1
pdf(j) = x*Ctq6Pdf(i,x,q)
enddo
else
mode = nset
call mrst2001(x,q,mode,upv,dnv,usea,dsea,str,chm,bot,glu)
pdf(-6) = 0
pdf(-5) = bot
pdf(-4) = chm
pdf(-3) = str
pdf(-2) = usea
pdf(-1) = dsea
pdf(0) = glu
pdf(1) = dnv + dsea
pdf(2) = upv + usea
pdf(3) = str
pdf(4) = chm
pdf(5) = bot
pdf(6) = 0
endif
pdf( 6) = 0
pdf(-6) = 0
return
end