hijing/src/pyhimult.f

0001
0002 C*********************************************************************
0003
0004       SUBROUTINE PYHIMULT(MMUL)
0005
0006 C...Initializes treatment of multiple interactions, selects kinematics
0007 C...of hardest interaction if low-pT physics included in run, and
0008 C...generates all non-hardest interactions.
0009       COMMON/LUJETS/N,K(9000,5),P(9000,5),V(9000,5)
0010       SAVE /LUJETS/
0011       COMMON/LUDAT1/MSTU(200),PARU(200),MSTJ(200),PARJ(200)
0012       SAVE /LUDAT1/
0013       COMMON/LUDAT2/KCHG(500,3),PMAS(500,4),PARF(2000),VCKM(4,4)
0014       SAVE /LUDAT2/
0015       COMMON/PYHISUBS/MSEL,MSUB(200),KFIN(2,-40:40),CKIN(200)
0016       SAVE /PYHISUBS/
0017       COMMON/PYHIPARS/MSTP(200),PARP(200),MSTI(200),PARI(200)
0018       SAVE /PYHIPARS/
0019       COMMON/PYHIINT1/MINT(400),VINT(400)
0020       SAVE /PYHIINT1/
0021       COMMON/PYHIINT2/ISET(200),KFPR(200,2),COEF(200,20),ICOL(40,4,2)
0022       SAVE /PYHIINT2/
0023       COMMON/PYHIINT3/XSFX(2,-40:40),ISIG(1000,3),SIGH(1000)
0024       SAVE /PYHIINT3/
0025       COMMON/PYHIINT5/NGEN(0:200,3),XSEC(0:200,3)
0026       SAVE /PYHIINT5/
0027       DIMENSION NMUL(20),SIGM(20),KSTR(500,2)
0028       SAVE XT2,XT2FAC,XC2,XTS,IRBIN,RBIN,NMUL,SIGM
0029
0030 C...Initialization of multiple interaction treatment.
0031       IF(MMUL.EQ.1) THEN
0032         IF(MSTP(122).GE.1) WRITE(MSTU(11),1000) MSTP(82)
0033         ISUB=96
0034         MINT(1)=96
0035         VINT(63)=0.
0036         VINT(64)=0.
0037         VINT(143)=1.
0038         VINT(144)=1.
0039
0040 C...Loop over phase space points: xT2 choice in 20 bins.
0041   100   SIGSUM=0.
0042         DO 120 IXT2=1,20
0043         NMUL(IXT2)=MSTP(83)
0044         SIGM(IXT2)=0.
0045         DO 110 ITRY=1,MSTP(83)
0046         RSCA=0.05*((21-IXT2)-RLU(0))
0047         XT2=VINT(149)*(1.+VINT(149))/(VINT(149)+RSCA)-VINT(149)
0048         XT2=MAX(0.01*VINT(149),XT2)
0049         VINT(25)=XT2
0050
0051 C...Choose tau and y*. Calculate cos(theta-hat).
0052         IF(RLU(0).LE.COEF(ISUB,1)) THEN
0053           TAUP=(2.*(1.+SQRT(1.-XT2))/XT2-1.)**RLU(0)
0054           TAU=XT2*(1.+TAUP)**2/(4.*TAUP)
0055         ELSE
0056           TAU=XT2*(1.+TAN(RLU(0)*ATAN(SQRT(1./XT2-1.)))**2)
0057         ENDIF
0058         VINT(21)=TAU
0059         CALL PYHIKLIM(2)
0060         RYST=RLU(0)
0061         MYST=1
0062         IF(RYST.GT.COEF(ISUB,7)) MYST=2
0063         IF(RYST.GT.COEF(ISUB,7)+COEF(ISUB,8)) MYST=3
0064         CALL PYHIKMAP(2,MYST,RLU(0))
0065         VINT(23)=SQRT(MAX(0.,1.-XT2/TAU))*(-1)**INT(1.5+RLU(0))
0066
0067 C...Calculate differential cross-section.
0068         VINT(71)=0.5*VINT(1)*SQRT(XT2)
0069         CALL PYHISIGH(NCHN,SIGS)
0070   110   SIGM(IXT2)=SIGM(IXT2)+SIGS
0071   120   SIGSUM=SIGSUM+SIGM(IXT2)
0072         SIGSUM=SIGSUM/(20.*MSTP(83))
0073
0074 C...Reject result if sigma(parton-parton) is smaller than hadronic one.
0075         IF(SIGSUM.LT.1.1*VINT(106)) THEN
0076           IF(MSTP(122).GE.1) WRITE(MSTU(11),1100) PARP(82),SIGSUM
0077           PARP(82)=0.9*PARP(82)
0078           VINT(149)=4.*PARP(82)**2/VINT(2)
0079           GOTO 100
0080         ENDIF
0081         IF(MSTP(122).GE.1) WRITE(MSTU(11),1200) PARP(82), SIGSUM
0082
0083 C...Start iteration to find k factor.
0084         YKE=SIGSUM/VINT(106)
0085         SO=0.5
0086         XI=0.
0087         YI=0.
0088         XK=0.5
0089         IIT=0
0090   130   IF(IIT.EQ.0) THEN
0091           XK=2.*XK
0092         ELSEIF(IIT.EQ.1) THEN
0093           XK=0.5*XK
0094         ELSE
0095           XK=XI+(YKE-YI)*(XF-XI)/(YF-YI)
0096         ENDIF
0097
0098 C...Evaluate overlap integrals.
0099         IF(MSTP(82).EQ.2) THEN
0100           SP=0.5*PARU(1)*(1.-EXP(-XK))
0101           SOP=SP/PARU(1)
0102         ELSE
0103           IF(MSTP(82).EQ.3) DELTAB=0.02
0104           IF(MSTP(82).EQ.4) DELTAB=MIN(0.01,0.05*PARP(84))
0105           SP=0.
0106           SOP=0.
0107           B=-0.5*DELTAB
0108   140     B=B+DELTAB
0109           IF(MSTP(82).EQ.3) THEN
0110             OV=EXP(-B**2)/PARU(2)
0111           ELSE
0112             CQ2=PARP(84)**2
0113             OV=((1.-PARP(83))**2*EXP(-MIN(100.,B**2))+2.*PARP(83)*
0114      &      (1.-PARP(83))*2./(1.+CQ2)*EXP(-MIN(100.,B**2*2./(1.+CQ2)))+
0115      &      PARP(83)**2/CQ2*EXP(-MIN(100.,B**2/CQ2)))/PARU(2)
0116           ENDIF
0117           PACC=1.-EXP(-MIN(100.,PARU(1)*XK*OV))
0118           SP=SP+PARU(2)*B*DELTAB*PACC
0119           SOP=SOP+PARU(2)*B*DELTAB*OV*PACC
0120           IF(B.LT.1..OR.B*PACC.GT.1E-6) GOTO 140
0121         ENDIF
0122         YK=PARU(1)*XK*SO/SP
0123
0124 C...Continue iteration until convergence.
0125         IF(YK.LT.YKE) THEN
0126           XI=XK
0127           YI=YK
0128           IF(IIT.EQ.1) IIT=2
0129         ELSE
0130           XF=XK
0131           YF=YK
0132           IF(IIT.EQ.0) IIT=1
0133         ENDIF
0134         IF(ABS(YK-YKE).GE.1E-5*YKE) GOTO 130
0135
0136 C...Store some results for subsequent use.
0137         VINT(145)=SIGSUM
0138         VINT(146)=SOP/SO
0139         VINT(147)=SOP/SP
0140
0141 C...Initialize iteration in xT2 for hardest interaction.
0142       ELSEIF(MMUL.EQ.2) THEN
0143         IF(MSTP(82).LE.0) THEN
0144         ELSEIF(MSTP(82).EQ.1) THEN
0145           XT2=1.
0146           XT2FAC=XSEC(96,1)/VINT(106)*VINT(149)/(1.-VINT(149))
0147         ELSEIF(MSTP(82).EQ.2) THEN
0148           XT2=1.
0149           XT2FAC=VINT(146)*XSEC(96,1)/VINT(106)*VINT(149)*(1.+VINT(149))
0150         ELSE
0151           XC2=4.*CKIN(3)**2/VINT(2)
0152           IF(CKIN(3).LE.CKIN(5).OR.MINT(82).GE.2) XC2=0.
0153         ENDIF
0154
0155       ELSEIF(MMUL.EQ.3) THEN
0156 C...Low-pT or multiple interactions (first semihard interaction):
0157 C...choose xT2 according to dpT2/pT2**2*exp(-(sigma above pT2)/norm)
0158 C...or (MSTP(82)>=2) dpT2/(pT2+pT0**2)**2*exp(-....).
0159         ISUB=MINT(1)
0160         IF(MSTP(82).LE.0) THEN
0161           XT2=0.
0162         ELSEIF(MSTP(82).EQ.1) THEN
0163           XT2=XT2FAC*XT2/(XT2FAC-XT2*LOG(RLU(0)))
0164         ELSEIF(MSTP(82).EQ.2) THEN
0165           IF(XT2.LT.1..AND.EXP(-XT2FAC*XT2/(VINT(149)*(XT2+
0166      &    VINT(149)))).GT.RLU(0)) XT2=1.
0167           IF(XT2.GE.1.) THEN
0168             XT2=(1.+VINT(149))*XT2FAC/(XT2FAC-(1.+VINT(149))*LOG(1.-
0169      &      RLU(0)*(1.-EXP(-XT2FAC/(VINT(149)*(1.+VINT(149)))))))-
0170      &      VINT(149)
0171           ELSE
0172             XT2=-XT2FAC/LOG(EXP(-XT2FAC/(XT2+VINT(149)))+RLU(0)*
0173      &      (EXP(-XT2FAC/VINT(149))-EXP(-XT2FAC/(XT2+VINT(149)))))-
0174      &      VINT(149)
0175           ENDIF
0176           XT2=MAX(0.01*VINT(149),XT2)
0177         ELSE
0178           XT2=(XC2+VINT(149))*(1.+VINT(149))/(1.+VINT(149)-
0179      &    RLU(0)*(1.-XC2))-VINT(149)
0180           XT2=MAX(0.01*VINT(149),XT2)
0181         ENDIF
0182         VINT(25)=XT2
0183
0184 C...Low-pT: choose xT2, tau, y* and cos(theta-hat) fixed.
0185         IF(MSTP(82).LE.1.AND.XT2.LT.VINT(149)) THEN
0186           IF(MINT(82).EQ.1) NGEN(0,1)=NGEN(0,1)-1
0187           IF(MINT(82).EQ.1) NGEN(ISUB,1)=NGEN(ISUB,1)-1
0188           ISUB=95
0189           MINT(1)=ISUB
0190           VINT(21)=0.01*VINT(149)
0191           VINT(22)=0.
0192           VINT(23)=0.
0193           VINT(25)=0.01*VINT(149)
0194
0195         ELSE
0196 C...Multiple interactions (first semihard interaction).
0197 C...Choose tau and y*. Calculate cos(theta-hat).
0198           IF(RLU(0).LE.COEF(ISUB,1)) THEN
0199             TAUP=(2.*(1.+SQRT(1.-XT2))/XT2-1.)**RLU(0)
0200             TAU=XT2*(1.+TAUP)**2/(4.*TAUP)
0201           ELSE
0202             TAU=XT2*(1.+TAN(RLU(0)*ATAN(SQRT(1./XT2-1.)))**2)
0203           ENDIF
0204           VINT(21)=TAU
0205           CALL PYHIKLIM(2)
0206           RYST=RLU(0)
0207           MYST=1
0208           IF(RYST.GT.COEF(ISUB,7)) MYST=2
0209           IF(RYST.GT.COEF(ISUB,7)+COEF(ISUB,8)) MYST=3
0210           CALL PYHIKMAP(2,MYST,RLU(0))
0211           VINT(23)=SQRT(MAX(0.,1.-XT2/TAU))*(-1)**INT(1.5+RLU(0))
0212         ENDIF
0213         VINT(71)=0.5*VINT(1)*SQRT(VINT(25))
0214
0215 C...Store results of cross-section calculation.
0216       ELSEIF(MMUL.EQ.4) THEN
0217         ISUB=MINT(1)
0218         XTS=VINT(25)
0219         IF(ISET(ISUB).EQ.1) XTS=VINT(21)
0220         IF(ISET(ISUB).EQ.2) XTS=(4.*VINT(48)+2.*VINT(63)+2.*VINT(64))/
0221      &  VINT(2)
0222         IF(ISET(ISUB).EQ.3.OR.ISET(ISUB).EQ.4) XTS=VINT(26)
0223         RBIN=MAX(0.000001,MIN(0.999999,XTS*(1.+VINT(149))/
0224      &  (XTS+VINT(149))))
0225         IRBIN=INT(1.+20.*RBIN)
0226         IF(ISUB.EQ.96) NMUL(IRBIN)=NMUL(IRBIN)+1
0227         IF(ISUB.EQ.96) SIGM(IRBIN)=SIGM(IRBIN)+VINT(153)
0228
0229 C...Choose impact parameter.
0230       ELSEIF(MMUL.EQ.5) THEN
0231         IF(MSTP(82).EQ.3) THEN
0232           VINT(148)=RLU(0)/(PARU(2)*VINT(147))
0233         ELSE
0234           RTYPE=RLU(0)
0235           CQ2=PARP(84)**2
0236           IF(RTYPE.LT.(1.-PARP(83))**2) THEN
0237             B2=-LOG(RLU(0))
0238           ELSEIF(RTYPE.LT.1.-PARP(83)**2) THEN
0239             B2=-0.5*(1.+CQ2)*LOG(RLU(0))
0240           ELSE
0241             B2=-CQ2*LOG(RLU(0))
0242           ENDIF
0243           VINT(148)=((1.-PARP(83))**2*EXP(-MIN(100.,B2))+2.*PARP(83)*
0244      &    (1.-PARP(83))*2./(1.+CQ2)*EXP(-MIN(100.,B2*2./(1.+CQ2)))+
0245      &    PARP(83)**2/CQ2*EXP(-MIN(100.,B2/CQ2)))/(PARU(2)*VINT(147))
0246         ENDIF
0247
0248 C...Multiple interactions (variable impact parameter) : reject with
0249 C...probability exp(-overlap*cross-section above pT/normalization).
0250         RNCOR=(IRBIN-20.*RBIN)*NMUL(IRBIN)
0251         SIGCOR=(IRBIN-20.*RBIN)*SIGM(IRBIN)
0252         DO 150 IBIN=IRBIN+1,20
0253         RNCOR=RNCOR+NMUL(IBIN)
0254   150   SIGCOR=SIGCOR+SIGM(IBIN)
0255         SIGABV=(SIGCOR/RNCOR)*VINT(149)*(1.-XTS)/(XTS+VINT(149))
0256         VINT(150)=EXP(-MIN(100.,VINT(146)*VINT(148)*SIGABV/VINT(106)))
0257
0258 C...Generate additional multiple semihard interactions.
0259       ELSEIF(MMUL.EQ.6) THEN
0260
0261 C...Reconstruct strings in hard scattering.
0262         ISUB=MINT(1)
0263         NMAX=MINT(84)+4
0264         IF(ISET(ISUB).EQ.1) NMAX=MINT(84)+2
0265         NSTR=0
0266         DO 170 I=MINT(84)+1,NMAX
0267         KCS=KCHG(LUCOMP(K(I,2)),2)*ISIGN(1,K(I,2))
0268         IF(KCS.EQ.0) GOTO 170
0269         DO 160 J=1,4
0270         IF(KCS.EQ.1.AND.(J.EQ.2.OR.J.EQ.4)) GOTO 160
0271         IF(KCS.EQ.-1.AND.(J.EQ.1.OR.J.EQ.3)) GOTO 160
0272         IF(J.LE.2) THEN
0273           IST=MOD(K(I,J+3)/MSTU(5),MSTU(5))
0274         ELSE
0275           IST=MOD(K(I,J+1),MSTU(5))
0276         ENDIF
0277         IF(IST.LT.MINT(84).OR.IST.GT.I) GOTO 160
0278         IF(KCHG(LUCOMP(K(IST,2)),2).EQ.0) GOTO 160
0279         NSTR=NSTR+1
0280         IF(J.EQ.1.OR.J.EQ.4) THEN
0281           KSTR(NSTR,1)=I
0282           KSTR(NSTR,2)=IST
0283         ELSE
0284           KSTR(NSTR,1)=IST
0285           KSTR(NSTR,2)=I
0286         ENDIF
0287   160   CONTINUE
0288   170   CONTINUE
0289
0290 C...Set up starting values for iteration in xT2.
0291         XT2=VINT(25)
0292         IF(ISET(ISUB).EQ.1) XT2=VINT(21)
0293         IF(ISET(ISUB).EQ.2) XT2=(4.*VINT(48)+2.*VINT(63)+2.*VINT(64))/
0294      &  VINT(2)
0295         IF(ISET(ISUB).EQ.3.OR.ISET(ISUB).EQ.4) XT2=VINT(26)
0296         ISUB=96
0297         MINT(1)=96
0298         IF(MSTP(82).LE.1) THEN
0299           XT2FAC=XSEC(ISUB,1)*VINT(149)/((1.-VINT(149))*VINT(106))
0300         ELSE
0301           XT2FAC=VINT(146)*VINT(148)*XSEC(ISUB,1)/VINT(106)*
0302      &    VINT(149)*(1.+VINT(149))
0303         ENDIF
0304         VINT(63)=0.
0305         VINT(64)=0.
0306         VINT(151)=0.
0307         VINT(152)=0.
0308         VINT(143)=1.-VINT(141)
0309         VINT(144)=1.-VINT(142)
0310
0311 C...Iterate downwards in xT2.
0312   180   IF(MSTP(82).LE.1) THEN
0313           XT2=XT2FAC*XT2/(XT2FAC-XT2*LOG(RLU(0)))
0314           IF(XT2.LT.VINT(149)) GOTO 220
0315         ELSE
0316           IF(XT2.LE.0.01*VINT(149)) GOTO 220
0317           XT2=XT2FAC*(XT2+VINT(149))/(XT2FAC-(XT2+VINT(149))*
0318      &    LOG(RLU(0)))-VINT(149)
0319           IF(XT2.LE.0.) GOTO 220
0320           XT2=MAX(0.01*VINT(149),XT2)
0321         ENDIF
0322         VINT(25)=XT2
0323
0324 C...Choose tau and y*. Calculate cos(theta-hat).
0325         IF(RLU(0).LE.COEF(ISUB,1)) THEN
0326           TAUP=(2.*(1.+SQRT(1.-XT2))/XT2-1.)**RLU(0)
0327           TAU=XT2*(1.+TAUP)**2/(4.*TAUP)
0328         ELSE
0329           TAU=XT2*(1.+TAN(RLU(0)*ATAN(SQRT(1./XT2-1.)))**2)
0330         ENDIF
0331         VINT(21)=TAU
0332         CALL PYHIKLIM(2)
0333         RYST=RLU(0)
0334         MYST=1
0335         IF(RYST.GT.COEF(ISUB,7)) MYST=2
0336         IF(RYST.GT.COEF(ISUB,7)+COEF(ISUB,8)) MYST=3
0337         CALL PYHIKMAP(2,MYST,RLU(0))
0338         VINT(23)=SQRT(MAX(0.,1.-XT2/TAU))*(-1)**INT(1.5+RLU(0))
0339
0340 C...Check that x not used up. Accept or reject kinematical variables.
0341         X1M=SQRT(TAU)*EXP(VINT(22))
0342         X2M=SQRT(TAU)*EXP(-VINT(22))
0343         IF(VINT(143)-X1M.LT.0.01.OR.VINT(144)-X2M.LT.0.01) GOTO 180
0344         VINT(71)=0.5*VINT(1)*SQRT(XT2)
0345         CALL PYHISIGH(NCHN,SIGS)
0346         IF(SIGS.LT.XSEC(ISUB,1)*RLU(0)) GOTO 180
0347
0348 C...Reset K, P and V vectors. Select some variables.
0349         DO 190 I=N+1,N+2
0350         DO 190 J=1,5
0351         K(I,J)=0
0352         P(I,J)=0.
0353   190   V(I,J)=0.
0354         RFLAV=RLU(0)
0355         PT=0.5*VINT(1)*SQRT(XT2)
0356         PHI=PARU(2)*RLU(0)
0357         CTH=VINT(23)
0358
0359 C...Add first parton to event record.
0360         K(N+1,1)=3
0361         K(N+1,2)=21
0362         IF(RFLAV.GE.MAX(PARP(85),PARP(86))) K(N+1,2)=
0363      &  1+INT((2.+PARJ(2))*RLU(0))
0364         P(N+1,1)=PT*COS(PHI)
0365         P(N+1,2)=PT*SIN(PHI)
0366         P(N+1,3)=0.25*VINT(1)*(VINT(41)*(1.+CTH)-VINT(42)*(1.-CTH))
0367         P(N+1,4)=0.25*VINT(1)*(VINT(41)*(1.+CTH)+VINT(42)*(1.-CTH))
0368         P(N+1,5)=0.
0369
0370 C...Add second parton to event record.
0371         K(N+2,1)=3
0372         K(N+2,2)=21
0373         IF(K(N+1,2).NE.21) K(N+2,2)=-K(N+1,2)
0374         P(N+2,1)=-P(N+1,1)
0375         P(N+2,2)=-P(N+1,2)
0376         P(N+2,3)=0.25*VINT(1)*(VINT(41)*(1.-CTH)-VINT(42)*(1.+CTH))
0377         P(N+2,4)=0.25*VINT(1)*(VINT(41)*(1.-CTH)+VINT(42)*(1.+CTH))
0378         P(N+2,5)=0.
0379
0380         IF(RFLAV.LT.PARP(85).AND.NSTR.GE.1) THEN
0381 C....Choose relevant string pieces to place gluons on.
0382           DO 210 I=N+1,N+2
0383           DMIN=1E8
0384           DO 200 ISTR=1,NSTR
0385           I1=KSTR(ISTR,1)
0386           I2=KSTR(ISTR,2)
0387           DIST=(P(I,4)*P(I1,4)-P(I,1)*P(I1,1)-P(I,2)*P(I1,2)-
0388      &    P(I,3)*P(I1,3))*(P(I,4)*P(I2,4)-P(I,1)*P(I2,1)-
0389      &    P(I,2)*P(I2,2)-P(I,3)*P(I2,3))/MAX(1.,P(I1,4)*P(I2,4)-
0390      &    P(I1,1)*P(I2,1)-P(I1,2)*P(I2,2)-P(I1,3)*P(I2,3))
0391           IF(ISTR.EQ.1.OR.DIST.LT.DMIN) THEN
0392             DMIN=DIST
0393             IST1=I1
0394             IST2=I2
0395             ISTM=ISTR
0396           ENDIF
0397   200     CONTINUE
0398
0399 C....Colour flow adjustments, new string pieces.
0400           IF(K(IST1,4)/MSTU(5).EQ.IST2) K(IST1,4)=MSTU(5)*I+
0401      &    MOD(K(IST1,4),MSTU(5))
0402           IF(MOD(K(IST1,5),MSTU(5)).EQ.IST2) K(IST1,5)=
0403      &    MSTU(5)*(K(IST1,5)/MSTU(5))+I
0404           K(I,5)=MSTU(5)*IST1
0405           K(I,4)=MSTU(5)*IST2
0406           IF(K(IST2,5)/MSTU(5).EQ.IST1) K(IST2,5)=MSTU(5)*I+
0407      &    MOD(K(IST2,5),MSTU(5))
0408           IF(MOD(K(IST2,4),MSTU(5)).EQ.IST1) K(IST2,4)=
0409      &    MSTU(5)*(K(IST2,4)/MSTU(5))+I
0410           KSTR(ISTM,2)=I
0411           KSTR(NSTR+1,1)=I
0412           KSTR(NSTR+1,2)=IST2
0413   210     NSTR=NSTR+1
0414
0415 C...String drawing and colour flow for gluon loop.
0416         ELSEIF(K(N+1,2).EQ.21) THEN
0417           K(N+1,4)=MSTU(5)*(N+2)
0418           K(N+1,5)=MSTU(5)*(N+2)
0419           K(N+2,4)=MSTU(5)*(N+1)
0420           K(N+2,5)=MSTU(5)*(N+1)
0421           KSTR(NSTR+1,1)=N+1
0422           KSTR(NSTR+1,2)=N+2
0423           KSTR(NSTR+2,1)=N+2
0424           KSTR(NSTR+2,2)=N+1
0425           NSTR=NSTR+2
0426
0427 C...String drawing and colour flow for q-qbar pair.
0428         ELSE
0429           K(N+1,4)=MSTU(5)*(N+2)
0430           K(N+2,5)=MSTU(5)*(N+1)
0431           KSTR(NSTR+1,1)=N+1
0432           KSTR(NSTR+1,2)=N+2
0433           NSTR=NSTR+1
0434         ENDIF
0435
0436 C...Update remaining energy; iterate.
0437         N=N+2
0438         IF(N.GT.MSTU(4)-MSTU(32)-10) THEN
0439           CALL LUERRM(11,'(PYHIMULT:) no more memory left in LUJETS')
0440           IF(MSTU(21).GE.1) RETURN
0441         ENDIF
0442         MINT(31)=MINT(31)+1
0443         VINT(151)=VINT(151)+VINT(41)
0444         VINT(152)=VINT(152)+VINT(42)
0445         VINT(143)=VINT(143)-VINT(41)
0446         VINT(144)=VINT(144)-VINT(42)
0447         IF(MINT(31).LT.240) GOTO 180
0448   220   CONTINUE
0449       ENDIF
0450
0451 C...Format statements for printout.
0452  1000 FORMAT(/1X,'****** PYHIMULT: initialization of multiple inter',
0453      &'actions for MSTP(82) =',I2,' ******')
0454  1100 FORMAT(8X,'pT0 =',F5.2,' GeV gives sigma(parton-parton) =',1P,
0455      &E9.2,' mb: rejected')
0456  1200 FORMAT(8X,'pT0 =',F5.2,' GeV gives sigma(parton-parton) =',1P,
0457      &E9.2,' mb: accepted')
0458
0459       RETURN
0460       END