pythiaeRHIC/pythia/pyinom.f

0001
0002 C*********************************************************************
0003
0004 C...PYINOM
0005 C...Finds the mass eigenstates and mixing matrices for neutralinos
0006 C...and charginos.
0007
0008       SUBROUTINE PYINOM
0009
0010 C...Double precision and integer declarations.
0011       IMPLICIT DOUBLE PRECISION(A-H, O-Z)
0012       IMPLICIT INTEGER(I-N)
0013       INTEGER PYCOMP
0014 C...Parameter statement to help give large particle numbers.
0015       PARAMETER (KSUSY1=1000000,KSUSY2=2000000,KTECHN=3000000,
0016      &KEXCIT=4000000,KDIMEN=5000000)
0017 C...Commonblocks.
0018       COMMON/PYDAT1/MSTU(200),PARU(200),MSTJ(200),PARJ(200)
0019       COMMON/PYDAT2/KCHG(500,4),PMAS(500,4),PARF(2000),VCKM(4,4)
0020       COMMON/PYMSSM/IMSS(0:99),RMSS(0:99)
0021       COMMON/PYSSMT/ZMIX(4,4),UMIX(2,2),VMIX(2,2),SMZ(4),SMW(2),
0022      &SFMIX(16,4),ZMIXI(4,4),UMIXI(2,2),VMIXI(2,2)
0023       SAVE /PYDAT1/,/PYDAT2/,/PYMSSM/,/PYSSMT/
0024
0025 C...Local variables.
0026       DOUBLE PRECISION XMW,XMZ,XM(4)
0027       DOUBLE PRECISION AR(4,4),WR(4),ZR(4,4),ZI(4,4),AI(4,4)
0028       DOUBLE PRECISION WI(4),FV1(4),FV2(4),FV3(4)
0029       DOUBLE PRECISION COSW,SINW
0030       DOUBLE PRECISION XMU
0031       DOUBLE PRECISION TANB,COSB,SINB
0032       DOUBLE PRECISION XM1,XM2,XM3,BETA
0033       DOUBLE PRECISION Q2,AEM,A1,A2,AQ,RM1,RM2
0034       DOUBLE PRECISION ARG,X0,X1,AX0,AX1,AT,BT
0035       DOUBLE PRECISION Y0,Y1,AMGX0,AM1X0,AMGX1,AM1X1
0036       DOUBLE PRECISION ARGX0,AR1X0,ARGX1,AR1X1
0037       DOUBLE PRECISION PYALPS,PYALEM
0038       DOUBLE PRECISION PYRNM3
0039       COMPLEX*16 CAR(4,4),CAI(4,4),CA1,CA2
0040       INTEGER IERR,INDEX(4),I,J,K,IOPT,ILR,KFNCHI(4)
0041       DATA KFNCHI/1000022,1000023,1000025,1000035/
0042
0043       IOPT=IMSS(2)
0044       IF(IMSS(1).EQ.2) THEN
0045         IOPT=1
0046       ENDIF
0047 C...M1, M2, AND M3 ARE INDEPENDENT
0048       IF(IOPT.EQ.0) THEN
0049         XM1=RMSS(1)
0050         XM2=RMSS(2)
0051         XM3=RMSS(3)
0052       ELSEIF(IOPT.GE.1) THEN
0053         Q2=PMAS(23,1)**2
0054         AEM=PYALEM(Q2)
0055         A2=AEM/PARU(102)
0056         A1=AEM/(1D0-PARU(102))
0057         XM1=RMSS(1)
0058         XM2=RMSS(2)
0059         IF(IMSS(1).EQ.2) XM1=RMSS(1)/RMSS(20)*A1*5D0/3D0
0060         IF(IOPT.EQ.1) THEN
0061           XM2=XM1*A2/A1*3D0/5D0
0062           RMSS(2)=XM2
0063         ELSEIF(IOPT.EQ.3) THEN
0064           XM1=XM2*5D0/3D0*A1/A2
0065           RMSS(1)=XM1
0066         ENDIF
0067         XM3=PYRNM3(XM2/A2)
0068         RMSS(3)=XM3
0069         IF(XM3.LE.0D0) THEN
0070           WRITE(MSTU(11),*) ' ERROR WITH M3 = ',XM3
0071           CALL PYSTOP(105)
0072         ENDIF
0073       ENDIF
0074
0075 C...GLUINO MASS
0076       IF(IMSS(3).EQ.1) THEN
0077         PMAS(PYCOMP(KSUSY1+21),1)=ABS(XM3)
0078       ELSE
0079         AQ=0D0
0080         DO 110 I=1,4
0081           DO 100 ILR=1,2
0082             RM1=PMAS(PYCOMP(ILR*KSUSY1+I),1)**2/XM3**2
0083             AQ=AQ+0.5D0*((2D0-RM1)*(RM1*LOG(RM1)-1D0)
0084      &      +(1D0-RM1)**2*LOG(ABS(1D0-RM1)))
0085   100     CONTINUE
0086   110   CONTINUE
0087
0088         DO 130 I=5,6
0089           DO 120 ILR=1,2
0090             RM1=PMAS(PYCOMP(ILR*KSUSY1+I),1)**2/XM3**2
0091             RM2=PMAS(I,1)**2/XM3**2
0092             ARG=(RM1-RM2-1D0)**2-4D0*RM2**2
0093             IF(ARG.GE.0D0) THEN
0094               X0=0.5D0*(1D0+RM2-RM1-SQRT(ARG))
0095               AX0=ABS(X0)
0096               X1=0.5D0*(1D0+RM2-RM1+SQRT(ARG))
0097               AX1=ABS(X1)
0098               IF(X0.EQ.1D0) THEN
0099                 AT=-1D0
0100                 BT=0.25D0
0101               ELSEIF(X0.EQ.0D0) THEN
0102                 AT=0D0
0103                 BT=-0.25D0
0104               ELSE
0105                 AT=0.5D0*LOG(ABS(1D0-X0))*(1D0-X0**2)+
0106      &          0.5D0*X0**2*LOG(AX0)
0107                 BT=(-1D0-2D0*X0)/4D0
0108               ENDIF
0109               IF(X1.EQ.1D0) THEN
0110                 AT=-1D0+AT
0111                 BT=0.25D0+BT
0112               ELSEIF(X1.EQ.0D0) THEN
0113                 AT=0D0+AT
0114                 BT=-0.25D0+BT
0115               ELSE
0116                 AT=0.5D0*LOG(ABS(1D0-X1))*(1D0-X1**2)+0.5D0*
0117      &          X1**2*LOG(AX1)+AT
0118                 BT=(-1D0-2D0*X1)/4D0+BT
0119               ENDIF
0120               AQ=AQ+AT+BT
0121             ELSE
0122               X0=0.5D0*(1D0+RM2-RM1)
0123               Y0=-0.5D0*SQRT(-ARG)
0124               AMGX0=SQRT(X0**2+Y0**2)
0125               AM1X0=SQRT((1D0-X0)**2+Y0**2)
0126               ARGX0=ATAN2(-X0,-Y0)
0127               AR1X0=ATAN2(1D0-X0,Y0)
0128               X1=X0
0129               Y1=-Y0
0130               AMGX1=AMGX0
0131               AM1X1=AM1X0
0132               ARGX1=ATAN2(-X1,-Y1)
0133               AR1X1=ATAN2(1D0-X1,Y1)
0134               AT=0.5D0*LOG(AM1X0)*(1D0-X0**2+3D0*Y0**2)
0135      &        +0.5D0*(X0**2-Y0**2)*LOG(AMGX0)
0136               BT=(-1D0-2D0*X0)/4D0+X0*Y0*( AR1X0-ARGX0 )
0137               AT=AT+0.5D0*LOG(AM1X1)*(1D0-X1**2+3D0*Y1**2)
0138      &        +0.5D0*(X1**2-Y1**2)*LOG(AMGX1)
0139               BT=BT+(-1D0-2D0*X1)/4D0+X1*Y1*( AR1X1-ARGX1 )
0140               AQ=AQ+AT+BT
0141             ENDIF
0142   120     CONTINUE
0143   130   CONTINUE
0144         PMAS(PYCOMP(KSUSY1+21),1)=ABS(XM3)*(1D0+PYALPS(XM3**2)
0145      &  /(2D0*PARU(2))*(15D0+AQ))
0146       ENDIF
0147
0148 C...NEUTRALINO MASSES
0149       DO 150 I=1,4
0150         DO 140 J=1,4
0151           AI(I,J)=0D0
0152   140   CONTINUE
0153   150 CONTINUE
0154       XMZ=PMAS(23,1)
0155       XMW=PMAS(24,1)
0156       XMU=RMSS(4)
0157       SINW=SQRT(PARU(102))
0158       COSW=SQRT(1D0-PARU(102))
0159       TANB=RMSS(5)
0160       BETA=ATAN(TANB)
0161       COSB=COS(BETA)
0162       SINB=TANB*COSB
0163
0164 C... Definitions:
0165 C...    psi^0 =(-i bino^0, -i wino^0, h_d^0(=H_1^0), h_u^0(=H_2^0))
0166 C... => L_neutralino = -1/2*(psi^0)^T * [AR] * psi^0 + h.c.
0167       AR(1,1) = XM1*COS(RMSS(30))
0168       AI(1,1) = XM1*SIN(RMSS(30))
0169       AR(2,2) = XM2*COS(RMSS(31))
0170       AI(2,2) = XM2*SIN(RMSS(31))
0171       AR(3,3) = 0D0
0172       AR(4,4) = 0D0
0173       AR(1,2) = 0D0
0174       AR(2,1) = 0D0
0175       AR(1,3) = -XMZ*SINW*COSB
0176       AR(3,1) = AR(1,3)
0177       AR(1,4) = XMZ*SINW*SINB
0178       AR(4,1) = AR(1,4)
0179       AR(2,3) = XMZ*COSW*COSB
0180       AR(3,2) = AR(2,3)
0181       AR(2,4) = -XMZ*COSW*SINB
0182       AR(4,2) = AR(2,4)
0183       AR(3,4) = -XMU*COS(RMSS(33))
0184       AI(3,4) = -XMU*SIN(RMSS(33))
0185       AR(4,3) = -XMU*COS(RMSS(33))
0186       AI(4,3) = -XMU*SIN(RMSS(33))
0187 C      CALL PYEIG4(AR,WR,ZR)
0188       CALL PYEICG(4,4,AR,AI,WR,WI,1,ZR,ZI,FV1,FV2,FV3,IERR)
0189       IF(IERR.NE.0) THEN
0190        WRITE(MSTU(11),*) ' PROBLEM WITH PYEICG IN PYINOM '
0191       ENDIF
0192       DO 160 I=1,4
0193         INDEX(I)=I
0194         XM(I)=ABS(WR(I))
0195   160 CONTINUE
0196       DO 180 I=2,4
0197         K=I
0198         DO 170 J=I-1,1,-1
0199           IF(XM(K).LT.XM(J)) THEN
0200             ITMP=INDEX(J)
0201             XTMP=XM(J)
0202             INDEX(J)=INDEX(K)
0203             XM(J)=XM(K)
0204             INDEX(K)=ITMP
0205             XM(K)=XTMP
0206             K=K-1
0207           ELSE
0208             GOTO 180
0209           ENDIF
0210   170   CONTINUE
0211   180 CONTINUE
0212
0213
0214       DO 210 I=1,4
0215         K=INDEX(I)
0216         SMZ(I)=WR(K)
0217         PMAS(PYCOMP(KFNCHI(I)),1)=ABS(SMZ(I))
0218         S=0D0
0219         DO 190 J=1,4
0220           S=S+ZR(J,K)**2+ZI(J,K)**2
0221   190   CONTINUE
0222         DO 200 J=1,4
0223           ZMIX(I,J)=ZR(J,K)/SQRT(S)
0224           ZMIXI(I,J)=ZI(J,K)/SQRT(S)
0225           IF(ABS(ZMIX(I,J)).LT.1D-6) ZMIX(I,J)=0D0
0226           IF(ABS(ZMIXI(I,J)).LT.1D-6) ZMIXI(I,J)=0D0
0227   200   CONTINUE
0228   210 CONTINUE
0229
0230 C...CHARGINO MASSES
0231 C.....Find eigenvectors of X X^*
0232       AI(1,1) = 0D0
0233       AI(2,2) = 0D0
0234       AR(1,1) = XM2**2+2D0*XMW**2*SINB**2
0235       AR(2,2) = XMU**2+2D0*XMW**2*COSB**2
0236       AR(1,2) = SQRT(2D0)*XMW*(XM2*COS(RMSS(31))*COSB+
0237      &XMU*COS(RMSS(33))*SINB)
0238       AI(1,2) = SQRT(2D0)*XMW*(XM2*SIN(RMSS(31))*COSB-
0239      &XMU*SIN(RMSS(33))*SINB)
0240       AR(2,1) = SQRT(2D0)*XMW*(XM2*COS(RMSS(31))*COSB+
0241      &XMU*COS(RMSS(33))*SINB)
0242       AI(2,1) = SQRT(2D0)*XMW*(-XM2*SIN(RMSS(31))*COSB+
0243      &XMU*SIN(RMSS(33))*SINB)
0244       CALL PYEICG(4,2,AR,AI,WR,WI,1,ZR,ZI,FV1,FV2,FV3,IERR)
0245       IF(IERR.NE.0) THEN
0246        WRITE(MSTU(11),*) ' PROBLEM WITH PYEICG IN PYINOM '
0247       ENDIF
0248       INDEX(1)=1
0249       INDEX(2)=2
0250       IF(WR(2).LT.WR(1)) THEN
0251         INDEX(1)=2
0252         INDEX(2)=1
0253       ENDIF
0254
0255       DO 240 I=1,2
0256         K=INDEX(I)
0257         SMW(I)=SQRT(WR(K))
0258         S=0D0
0259         DO 220 J=1,2
0260           S=S+ZR(J,K)**2+ZI(J,K)**2
0261   220   CONTINUE
0262         DO 230 J=1,2
0263           UMIX(I,J)=ZR(J,K)/SQRT(S)
0264           UMIXI(I,J)=-ZI(J,K)/SQRT(S)
0265           IF(ABS(UMIX(I,J)).LT.1D-6) UMIX(I,J)=0D0
0266           IF(ABS(UMIXI(I,J)).LT.1D-6) UMIXI(I,J)=0D0
0267   230   CONTINUE
0268   240 CONTINUE
0269 C...Force chargino mass > neutralino mass
0270       IF(ABS(SMW(1)).LT.ABS(SMZ(1))+2D0*PMAS(PYCOMP(111),1)) THEN
0271         CALL PYERRM(18,'(PYINOM:) '//
0272      &      'forcing m(~chi+_1) > m(~chi0_1) + 2m(pi0)')
0273         SMW(1)=SIGN(ABS(SMZ(1))+2D0*PMAS(PYCOMP(111),1),SMW(1))
0274       ENDIF
0275       PMAS(PYCOMP(KSUSY1+24),1)=SMW(1)
0276       PMAS(PYCOMP(KSUSY1+37),1)=SMW(2)
0277
0278 C.....Find eigenvectors of X^* X
0279       AI(1,1) = 0D0
0280       AI(2,2) = 0D0
0281       AR(1,1) = XM2**2+2D0*XMW**2*COSB**2
0282       AR(2,2) = XMU**2+2D0*XMW**2*SINB**2
0283       AR(1,2) = SQRT(2D0)*XMW*(XM2*COS(RMSS(31))*SINB+
0284      &XMU*COS(RMSS(33))*COSB)
0285       AI(1,2) = SQRT(2D0)*XMW*(-XM2*SIN(RMSS(31))*SINB+
0286      &XMU*SIN(RMSS(33))*COSB)
0287       AR(2,1) = SQRT(2D0)*XMW*(XM2*COS(RMSS(31))*SINB+
0288      &XMU*COS(RMSS(33))*COSB)
0289       AI(2,1) = SQRT(2D0)*XMW*(XM2*SIN(RMSS(31))*SINB-
0290      &XMU*SIN(RMSS(33))*COSB)
0291       CALL PYEICG(4,2,AR,AI,WR,WI,1,ZR,ZI,FV1,FV2,FV3,IERR)
0292       IF(IERR.NE.0) THEN
0293        WRITE(MSTU(11),*) ' PROBLEM WITH PYEICG IN PYINOM '
0294       ENDIF
0295       INDEX(1)=1
0296       INDEX(2)=2
0297       IF(WR(2).LT.WR(1)) THEN
0298         INDEX(1)=2
0299         INDEX(2)=1
0300       ENDIF
0301
0302       DO 270 I=1,2
0303         K=INDEX(I)
0304         S=0D0
0305         DO 250 J=1,2
0306           S=S+ZR(J,K)**2+ZI(J,K)**2
0307   250   CONTINUE
0308         DO 260 J=1,2
0309           VMIX(I,J)=ZR(J,K)/SQRT(S)
0310           VMIXI(I,J)=-ZI(J,K)/SQRT(S)
0311           IF(ABS(VMIX(I,J)).LT.1D-6) VMIX(I,J)=0D0
0312           IF(ABS(VMIXI(I,J)).LT.1D-6) VMIXI(I,J)=0D0
0313   260   CONTINUE
0314   270 CONTINUE
0315
0316
0317       RETURN
0318       END