SUBROUTINE ESC96(R,TYPE,EFAC,FI0,FI1,FI2,VPOT) ************************************************************************ ** Nijmegen Extended Soft-Core (ESC) NN potential in coordinate space ** including pi-meson BW, TMO, and pair potentials ** Version 1.1: October 1999 ** Email: info@nn-online.org ** WWW: http://nn-online.org ** References: Rijken & Stoks, Phys. Rev. C 54 (1996), 2851, 2869. ** ** INPUT: R in fermi ** ----- TYPE 'PP', 'NP', 'PN', or 'NN' (character*2) ** EFAC M/E factor for energy dependent one-pion exchange ** ** OUTPUT: VPOT 2x2 potential matrix in MeV on LSJ-basis ** ------ FI0,FI1,FI2 2x2 nonlocal contribution and its derivatives ** with dimensions [0], [MeV], [MeV**2], ** including the tensor nonlocal contribution. ** ** COMMON/EMANHP/PHNAME has to be filled beforehand. ** PHNAME is character*3 and contains the name of the phase ** shift in the spectral notation. ** - singlets: 1S0 1P1 1D2 1F3 1G4 ... ** - triplets uncoupled: 3P0 3P1 3D2 3F3 3G4 ... ** - triplets coupled: 3C1 3C2 3C3 3C4 ... ** where 3C1 denotes 3S1 - epsilon1 - 3D1 channel ** 3C2 denotes 3P2 - epsilon2 - 3F2 channel ... ** DATA statements: ** IPAIR = 1: pair potentials included ** IPIMES= 1: pi-meson potentials included ** IAXIAL= 1: axial OBE potentials included ** NLOC = 1: standard nonlocalities ** 2: nonlocalities in spin-spin and tensor of OPE ************************************************************************ IMPLICIT REAL*8 (A-H,O-Z) CHARACTER PHNAME*3,PHNAM0*3, TYPE*2,TYP0*2 LOGICAL FIRST, CALC, NWGRID, niso0,niso1 INTEGER SPIN REAL*8 VPOT(2,2), FI0(2,2),FI1(2,2),FI2(2,2) PARAMETER(IPTMAX=2000) COMMON/EMANHP/PHNAME COMMON/POT2/ UC,USS,UT,ULS,ULSA,UQ12,FIC0,FIC1,FIC2 COMMON/VOPE/ VPIS,VPIT,VPIQ,FPIS0,FPIS1,FPIS2,FPIT0,FPIT1,FPIT2 COMMON/VENER/ PLMEV2,VEPIS,VEPIT COMMON/VPIMES/VCPI,VSPI,VTPI,VLPI, VCPM,VSPM,VTPM,VLPM COMMON/VPAIR/ V1C,V1S,V1T,V1L, V2C,V2S,V2T,V2L COMMON/POTFSV/ICAL,IPT,POTMAT(IPTMAX,10),OPEMAT(IPTMAX,2) DATA NLOC/1/, NADIA/0/, IPAIR/1/, IPIMES/1/, IAXIAL/0/ DATA NOWRIT/1/ DATA R0/-1D0/, ISO0/-1/, PHNAM0/'***'/, TYP0/'XX'/ data niso0/.true./, niso1/.true./, rw0/1d-3/, eps/1d-4/ SAVE FIRST,CALC,NWGRID, NCHAN,SPIN,L,J,ISO, TJ,TJM,TJP,SQJJ DO 50 I=1,2*2 VPOT(I,1)=0D0 FI0(I,1)=0D0 FI1(I,1)=0D0 FI2(I,1)=0D0 50 CONTINUE IF(PHNAME.NE.PHNAM0) THEN IF(NOWRIT.EQ.1) PHNAM0=PHNAME NCHAN=1 IF(PHNAME(2:2).EQ.'C') NCHAN=2 IF(PHNAME(1:1).EQ.'1') SPIN=0 IF(PHNAME(1:1).EQ.'3') SPIN=1 READ(PHNAME,'(2X,I1)') J L=J IF(PHNAME.EQ.'3P0') L=1 IF(NCHAN.EQ.2) L=J-1 ISO=MOD(SPIN+L+1,2) TJ=2D0*J+1D0 TJM=(J-1D0) TJP=(J+2D0) SQJJ=DSQRT(J*(J+1D0)) ENDIF IF(NOWRIT.EQ.1) THEN *** Potential is not stored in POTMAT CALL NYMISO(R,TYPE,ISO,NLOC,NADIA) VC =UC VSS =USS VT =UT VLS =ULS VLSA=ULSA VQ12=UQ12 IF(IPAIR.EQ.1) THEN CALL PAIR(R,TYPE,ISO) VC =VC + (V1C+V2C) VSS=VSS+ (V1S+V2S) VT =VT + (V1T+V2T) VLS=VLS+ (V1L+V2L) ENDIF IF(IPIMES.EQ.1) THEN CALL PIMES(R,TYPE,ISO) VC =VC + (VCPI+VCPM) VSS=VSS+ (VSPI+VSPM) VT =VT + (VTPI+VTPM) VLS=VLS+ (VLPI+VLPM) ENDIF IF(IAXIAL.EQ.1) THEN CALL AXIAL(R,ISO,SPIN,VAXSS,VAXT,VAXLS,FI0AX,FI1AX,FI2AX) VSS =VSS +VAXSS VT =VT +VAXT VLS =VLS +VAXLS FIC0=FIC0+FI0AX FIC1=FIC1+FI1AX FIC2=FIC2+FI2AX ENDIF GOTO 13 ENDIF NWGRID=.FALSE. IF(PHNAME.NE.PHNAM0 .AND. . (PHNAME.EQ.'1S0' .OR. PHNAM0.EQ.'1S0')) NWGRID=.TRUE. IF(ISO.NE.ISO0 .OR. TYPE.NE.TYP0 .OR. NWGRID) FIRST=.TRUE. IF(ICAL.EQ.0) FIRST=.TRUE. IF(FIRST) THEN CALC=.TRUE. IPT=1 ENDIF IF(.NOT.FIRST .AND. R.LT.R0) THEN CALC=.FALSE. IPT=1 ENDIF ICAL=1 FIRST=.FALSE. ISO0=ISO TYP0=TYPE PHNAM0=PHNAME R0=R IF(.NOT.CALC) THEN 11 R0=POTMAT(IPT,1) IF(R0.LT.R) THEN IPT=IPT+1 IF(IPT.GT.IPTMAX) GOTO 12 GOTO 11 ELSEIF(DABS(R0-R).LT.1D-6) THEN VC =POTMAT(IPT,2) VSS =POTMAT(IPT,3) VT =POTMAT(IPT,4) VLS =POTMAT(IPT,5) VLSA=POTMAT(IPT,6) VQ12=POTMAT(IPT,7) FIC0=POTMAT(IPT,8) FIC1=POTMAT(IPT,9) FIC2=POTMAT(IPT,10) c FIS0=POTMAT(IPT,11) c FIS1=POTMAT(IPT,12) c FIS2=POTMAT(IPT,13) c FIT0=POTMAT(IPT,14) c FIT1=POTMAT(IPT,15) c FIT2=POTMAT(IPT,16) c IF(NADIA.EQ.1) VEPIS=POTMAT(IPT,17) c IF(NADIA.EQ.1) VEPIT=POTMAT(IPT,18) VPIS=OPEMAT(IPT,1) VPIT=OPEMAT(IPT,2) GOTO 13 ENDIF 12 STOP '*** NIJMESC96: stored grid file differs with input' ELSE CALL NYMISO(R,TYPE,ISO,NLOC,NADIA) VC =UC VSS =USS VT =UT VLS =ULS VLSA=ULSA VQ12=UQ12 IF(IPAIR.EQ.1) THEN CALL PAIR(R,TYPE,ISO) VC =VC + (V1C+V2C) VSS=VSS+ (V1S+V2S) VT =VT + (V1T+V2T) VLS=VLS+ (V1L+V2L) ENDIF IF(IPIMES.EQ.1) THEN CALL PIMES(R,TYPE,ISO) VC =VC + (VCPI+VCPM) VSS=VSS+ (VSPI+VSPM) VT =VT + (VTPI+VTPM) VLS=VLS+ (VLPI+VLPM) ENDIF IF(IAXIAL.EQ.1) THEN CALL AXIAL(R,ISO,SPIN,VAXSS,VAXT,VAXLS,FI0AX,FI1AX,FI2AX) VSS =VSS +VAXSS VT =VT +VAXT VLS =VLS +VAXLS FIC0=FIC0+FI0AX FIC1=FIC1+FI1AX FIC2=FIC2+FI2AX ENDIF POTMAT(IPT,1) =R POTMAT(IPT,2) =VC POTMAT(IPT,3) =VSS POTMAT(IPT,4) =VT POTMAT(IPT,5) =VLS POTMAT(IPT,6) =VLSA POTMAT(IPT,7) =VQ12 POTMAT(IPT,8) =FIC0 POTMAT(IPT,9) =FIC1 POTMAT(IPT,10)=FIC2 c POTMAT(IPT,11)=FIS0 c POTMAT(IPT,12)=FIS1 c POTMAT(IPT,13)=FIS2 c POTMAT(IPT,14)=FIT0 c POTMAT(IPT,15)=FIT1 c POTMAT(IPT,16)=FIT2 c IF(NADIA.EQ.1) POTMAT(IPT,17)=VEPIS c IF(NADIA.EQ.1) POTMAT(IPT,18)=VEPIT OPEMAT(IPT,1)=VPIS OPEMAT(IPT,2)=VPIT IPT=IPT+1 IF(IPT.GT.IPTMAX) STOP '*** NIJMESC96: POTMAT is too small' ENDIF 13 CONTINUE IF(NADIA.EQ.1) THEN *** Non-adiabatic contribution to OPE VSS=VSS+VPIS + PLMEV2*VEPIS VT =VT +VPIT + PLMEV2*VEPIT ELSEIF(NLOC.EQ.1) THEN *** Energy-dependent (M/E-factor) one-poin exchange c CALL PIOME(R,TYPE,EFAC,VPIS,VPIT) VSS=VSS+VPIS*EFAC VT =VT +VPIT*EFAC ENDIF IF(NCHAN.EQ.1) THEN IF(SPIN.EQ.0) THEN VPOT(1,1) = VC - 3D0*VSS - J*(J+1)*VQ12 ELSEIF(L.EQ.J) THEN VPOT(1,1) = VC + VSS + 2D0*VT - VLS + (1D0-J*(J+1))*VQ12 ELSEIF(PHNAME.EQ.'3P0') THEN VPOT(1,1) = VC + VSS - 4D0*VT - 2D0*VLS + 4D0*VQ12 ENDIF ELSE VPOT(1,1) = VC + VSS - 2D0*(J-1)/(2*J+1)*VT + . (J-1)*VLS + (J-1)*(J-1)*VQ12 VPOT(2,2) = VC + VSS - 2D0*(J+2)/(2*J+1)*VT - . (J+2)*VLS + (J+2)*(J+2)*VQ12 VPOT(1,2) = 6D0*DSQRT(J*(J+1)*1D0)/(2*J+1)*VT VPOT(2,1) = VPOT(1,2) ENDIF IF(NLOC.EQ.1 .AND. NADIA.EQ.0) THEN FI0(1,1) = FIC0 FI1(1,1) = FIC1 FI2(1,1) = FIC2 IF(NCHAN.EQ.2) THEN FI0(2,2) = FIC0 FI1(2,2) = FIC1 FI2(2,2) = FIC2 ENDIF c ELSEIF(NLOC.EQ.2 .OR. NADIA.EQ.1) THEN c IF(NCHAN.EQ.1) THEN c IF(SPIN.EQ.0) THEN c FI0(1,1) = FIC0 - 3D0*FIS0 c FI1(1,1) = FIC1 - 3D0*FIS1 c FI2(1,1) = FIC2 - 3D0*FIS2 c ELSEIF(L.EQ.J) THEN c FI0(1,1) = FIC0 + FIS0 + 2D0*FIT0 c FI1(1,1) = FIC1 + FIS1 + 2D0*FIT1 c FI2(1,1) = FIC2 + FIS2 + 2D0*FIT2 c ELSEIF(PHNAME.EQ.'3P0') THEN c FI0(1,1) = FIC0 + FIS0 - 4D0*FIT0 c FI1(1,1) = FIC1 + FIS1 - 4D0*FIT1 c FI2(1,1) = FIC2 + FIS2 - 4D0*FIT2 c ENDIF c ELSEIF(NCHAN.EQ.2) THEN c FI0(1,1) = FIC0 + FIS0 - 2D0*TJM/TJ*FIT0 c FI1(1,1) = FIC1 + FIS1 - 2D0*TJM/TJ*FIT1 c FI2(1,1) = FIC2 + FIS2 - 2D0*TJM/TJ*FIT2 c FI0(1,2) = 6D0*SQJJ/TJ*FIT0 c FI1(1,2) = 6D0*SQJJ/TJ*FIT1 c FI2(1,2) = 6D0*SQJJ/TJ*FIT2 c FI0(2,1) = 6D0*SQJJ/TJ*FIT0 c FI1(2,1) = 6D0*SQJJ/TJ*FIT1 c FI2(2,1) = 6D0*SQJJ/TJ*FIT2 c FI0(2,2) = FIC0 + FIS0 - 2D0*TJP/TJ*FIT0 c FI1(2,2) = FIC1 + FIS1 - 2D0*TJP/TJ*FIT1 c FI2(2,2) = FIC2 + FIS2 - 2D0*TJP/TJ*FIT2 c ENDIF ENDIF if( ((niso0.and.iso.eq.0) .or. (niso1.and.iso.eq.1)) . .and. dabs(r-rw0).lt.eps) then eptf=1d0+2d0*fi0(1,1) v11=vc-3d0*vss-j*(j+1)*vq12 v33=vc+vss+2d0*vt-vls+(1d0-j*(j+1))*vq12 write(*,'(a6,2x,d15.8,2(a14,f12.5))') 'r(fm)=',r,'eff-sing=', . v11/eptf-(fi1(1,1)/eptf)**2/939d0,'eff-trip=', . v33/eptf-(fi1(1,1)/eptf)**2/939d0 write(*,'(a4,i1,5x,5f13.5)') 'iso=',iso, . vc,vss,vt,vls,vq12 write(*,'(6x,a3,1x,5f13.5)') 'OBE', . uc,uss+vpis*efac,ut+vpit*efac,uls,uq12+vpiq*efac write(*,99) ' pi-pi',vcpi,vspi,vtpi,vlpi write(*,99) 'pi-mes',vcpm,vspm,vtpm,vlpm write(*,99) ' pair',v1c+v2c,v1s+v2s,v1t+v2t,v1l+v2l write(*,99) '1-pair',v1c,v1s,v1t,v1l write(*,99) '2-pair',v2c,v2s,v2t,v2l if(iso.eq.0) niso0=.false. if(iso.eq.1) niso1=.false. 99 format(3x,a6,1x,4f13.5) endif RETURN END ************************************************************************ ** NIJMEGEN NUCLEON - NUCLEON POTENTIAL PROGRAM ** ** ----------------------------------------------- ** ************************************************************************ SUBROUTINE NYMISO(R,TYPE,ISO,NLOC,NADIA) IMPLICIT REAL*8 (A-H,O-Z) REAL*8 NEUTM CHARACTER PHNAME*3,PHNAM0*3, TYPE*2,TYP0*2 LOGICAL FIRST, BRDEP,BRDRO,BRDROC common/hoiho/ scal COMMON/EMANHP/PHNAME COMMON/MESONM/AMPI,AMETA,AMETP, AMRO,AMOM,AMFI, AMA0,AMEP,AMF0, . AMPIC,AMROC,AMA0C,AWPIC,AWROC,AWA0C, AVSC COMMON/MEANM/ PIM,ROM,A0M COMMON/PARAMS/PAR(8,5) COMMON/COPLNG/ALPV,THPV,PV1, FPI,FETA,FETP,FPI2,FETA2,FETP2, . ALVE,THV ,GV1, GRO,GOM,GFI, GRO2,GOM2,GFI2, . ALVM, FV1, FRO,FOM,FFI, FRO2,FOM2,FFI2, . ALGS,THGS,GS1, GA0,GEP,GF0, GA02,GEP2,GF02, . GFPRO,GFPOM,GFPFI, GFMRO,GFMOM,GFMFI, . FPIC2,GROC2,FROC2,GFPROC,GFMROC,GA0C2 COMMON/KAPPAV/AKRO,AKOM,AKFI COMMON/YUKEFF/ARO,AMR1,AROC,AMRC1,AWRC1,BRO,AMR2,BROC,AMRC2,AWRC2, . AVSC1,AVSC2, AEPS,AME1,BEPS,AME2 COMMON/BRDRHO/GAMRO,THRRO,NRO,BRDRO,GAMROC,THRROC,NROC,BRDROC COMMON/BRDEPS/GAMEP,THREP0,THREPC,NEP,BRDEP COMMON/SCALIN/AMT,AMPV,AMT24 COMMON/AMSCAL/PIOMS2,PIOMS COMMON/CUTOFF/ALAM1,ALAMP1,ALAMV1,ALAMS1, . ALAM0,ALAMP0,ALAMV0,ALAMS0 COMMON/CUTTEN/ALAMT1,ALAVT1,ALAMT0,ALAVT0 COMMON/POMRON/GPOM,FPOM2,AMPOM,AMPOM2,AMPOM4, . GA2D,FA2D2,AMA2D,AMA2D2,AMA2D4 COMMON/AMCOEF/ALF,REDM, AMY,AMY2,AMYI,AMYI2, AMN,AMN2,AMNI,AMNI2, . AMYPN,AMYMN, AMYPN2,AMYMN2, AY2PN2,AY2MN2, . AMYN,AMYNI,AMYNI2, AYPNI2,AYMNI2 COMMON/GPAIR/ GPIETA,GPIETP,GSISI, GPISIG,GPIPOM, . GPIPI0,GPIPI1,FPIPI1,GPIRO0,GPIRO1,GPIOM COMMON/HBARC/ HBC DATA FPPPI0/0.0745D0/, FIRST/.TRUE./ DATA TYP0/'XX'/, PHNAM0/'***'/, RGT/25D0/ DATA PI/3.14159265358979323846D0/, SRPI/1.7724538509055160273D0/ SAVE ISIGN,CONV,PROTM,NEUTM IF(.NOT.FIRST) GOTO 10 FIRST = .FALSE. CONV = PI/180D0 ** Nucleon and meson masses (Particle Data Group 1992) ** Only include meson-mass difference for pions HBC = 197.327053D0 PROTM = 938.27231D0 NEUTM = 939.56563D0 AMT = 938.27231D0 AMT24 = 4D0*AMT*AMT AMPV = 139.5679D0 AMPIC = 139.5679D0 AMPI = 134.9743D0 PIM = (AMPI+2D0*AMPIC)/3D0 AMETA = 547.45D0 AMETP = 957.75D0 ROM = 768.1D0 AMROC = ROM AMRO = ROM AMOM = 781.95D0 AMFI =1019.413D0 A0M = 982.7D0 AMA0C = A0M AMA0 = A0M AMEP = 760D0 AMF0 = 974.1D0 AVSC = 1D0/AMROC**2 FAC = 1D0 AWPIC = AMPIC AWROC = AMROC AWA0C = AMA0C PIOMS = AMPV PIOMS2= AMPV*AMPV ** Broad rho meson: spectral density to two effective Yukawa's ** Yukawa's fitted to PHI0C 0.001 - 2 fm (steps 0.005, ALAM1=850 MeV) BRDRO = .FALSE. GAMRO = 151.5D0 THRRO = PIM+PIM NRO = 1 ARO = 0.35765D0 AMR1 = 667.877D0 BRO = 0.47396D0 AMR2 = 914.234D0 BRDROC= BRDRO GAMROC= GAMRO THRROC= THRRO NROC = NRO AROC = ARO AMRC1 = AMR1 AWRC1 = AMRC1 BROC = BRO AMRC2 = AMR2 AWRC2 = AMRC2 AVSC1 = 1D0/AMRC1**2 AVSC2 = 1D0/AMRC2**2 ** Broad epsilon meson: spectral density to two effective Yukawa's ** Yukawa's fitted to PHI0C 0.001 - 2 fm (steps 0.005, ALAM0=800 MeV) BRDEP = .FALSE. GAMEP = 800D0 THREP0= AMPI+AMPI THREPC= AMPIC+AMPIC NEP = 0 AEPS = 0.16914D0 AME1 = 487.913D0 BEPS = 0.61302D0 AME2 =1021.402D0 CALL ESCPAR C pseudovec vector tensor scalar pomeron C FPI GRO AKRO GA0 GA2D C PV1 GV1 FV1 GEP GPOM C ALPV ALVE ALVM GF0 AMPOM C THPV THV GA0C% THGS AMA2D C Cut-off ALAMP1 ALAMV1 ALAM1 ALAMS1 SCAL C ALAMP0 ALAMV0 ALAM0 ALAMS0 GPIOM C Pairs GPIETA GPIETP GSISI GPISIG GPIPOM C GPIPI0 GPIPI1 FPIPI1 GPIRO0 GPIRO1 ** (joined) cutoffs for pseudoscalar, vector, scalar ALAMP1 = PAR(5,1) ALAMV1 = PAR(5,2) ALAMS1 = PAR(5,4) ALAM1 = PAR(5,3) ALAMP0 = PAR(6,1) ALAMV0 = PAR(6,2) ALAMS0 = PAR(6,4) ALAM0 = PAR(6,3) IF(ALAM1.NE.0D0) THEN ALAMP1=ALAM1 ALAMV1=ALAM1 ALAMS1=ALAM1 IF(ALAM0.EQ.0D0) ALAM0=ALAM1 ENDIF IF(ALAM0.NE.0D0) THEN ALAMP0=ALAM0 ALAMV0=ALAM0 ALAMS0=ALAM0 ENDIF alamt1=alamv1 alavt1=alamv1 alamt0=alamv0 alavt0=alamv0 ** pseudovector couplings c FPI = PAR(1,1) FPI = DSQRT(FPPPI0*DEXP(-(AMPI/ALAMP1)**2)) PAR(1,1)=FPI PV1 = PAR(2,1) ALPV = PAR(3,1) THPV = PAR(4,1)*CONV ** vector couplings GRO = PAR(1,2) GV1 = PAR(2,2) ALVE = PAR(3,2) THV = PAR(4,2)*CONV ** tensor couplings AKRO = PAR(1,3) FV1 = PAR(2,3) ALVM = PAR(3,3) ** scalar couplings GA0 = PAR(1,4) GEP = PAR(2,4) GF0 = PAR(3,4) THGS = PAR(4,4)*CONV ** diffractive contribution GA2D = PAR(1,5) GPOM = PAR(2,5) AMPOM = PAR(3,5) AMA2D = PAR(4,5) IF(AMA2D.EQ.0D0) AMA2D=AMPOM AMPOM2= AMPOM*AMPOM AMPOM4= AMPOM2*AMPOM2 AMA2D2= AMA2D*AMA2D AMA2D4= AMA2D2*AMA2D2 ** useful parameter scal = PAR(5,5) ** Pair couplings GPIETA = PAR(7,1) GPIETP = PAR(7,2) GSISI = PAR(7,3) GPISIG = PAR(7,4) GPIPOM = PAR(7,5) GPIPI0 = PAR(8,1) GPIPI1 = PAR(8,2) FPIPI1 = PAR(8,3) GPIRO0 = PAR(8,4) GPIRO1 = PAR(8,5) GPIOM = PAR(6,5) ** define the physical coupling constants CALL NYMCPFF CALL CHIRAL 10 CONTINUE IF(TYPE.EQ.TYP0) GOTO 20 TYP0=TYPE IF(TYPE.EQ.'PP') THEN AMY = PROTM AMN = PROTM I3Y = 1 I3N = 1 ELSEIF(TYPE.EQ.'NP') THEN AMY = NEUTM AMN = PROTM I3Y =-1 I3N = 1 ELSEIF(TYPE.EQ.'PN') THEN AMY = PROTM AMN = NEUTM I3Y = 1 I3N =-1 ELSEIF(TYPE.EQ.'NN') THEN AMY = NEUTM AMN = NEUTM I3Y =-1 I3N =-1 ENDIF CALL AMFACNN ISIGN= I3Y*I3N FPI2 = DABS(FPI2)*ISIGN GRO2 = DABS(GRO2)*ISIGN FRO2 = DABS(FRO2)*ISIGN GA02 = DABS(GA02)*ISIGN GFPRO= DABS(GFPRO)*ISIGN GFMRO= DABS(GFMRO)*ISIGN ** isospin dependence 20 FACISO= 4D0*ISO-2D0 GA0C2 = FACISO*GA0*GA0 FPIC2 = FACISO*FPI*FPI/(AMPV*AMPV) GROC2 = FACISO*GRO*GRO FROC2 = FACISO*FRO*FRO/AMT24 GFPROC= FACISO*(GRO*FRO+FRO*GRO)/(2D0*AMT) GFMROC= FACISO*(GRO*FRO-FRO*GRO)/(2D0*AMT) ** Phenomenological phase-shift difference in PP and NP 1S0 IF(PHNAME.EQ.'1S0' .AND. (TYPE.EQ.'NP'.OR.TYPE.EQ.'PN')) . GA0C2=GA0C2*(1D0+PAR(4,3)) FPOM2 = (GPOM*AMPOM/AMT)**2 * 4D0*AMPOM/SRPI FA2D2 = (GA2D*AMA2D/AMT)**2 * 4D0*AMA2D/SRPI RH = R/HBC RGT= R CALL RSPACFF(RH,TYPE,NADIA) CALL POTFNN(RH,TYPE,ISIGN,FACISO,NLOC,NADIA) RETURN END ************************************************************************ SUBROUTINE ESCPAR IMPLICIT REAL*8 (A-H,O-Z) COMMON/PARAMS/PAR(8,5) DIMENSION PARTR(5,8) DATA PARTR/ . .2695464D+00,.6198737D+00,.3705890D+01,.1868811D+01,.0000000D+00 .,.1632696D+00,.3419249D+01,-.2432022D+00,.5810781D+00,.1954388D+01 .,.3550000D+00,.1000000D+01,.400000D+00,-.2457693D+01,.3091017D+03 .,-.230000D+02,.3750000D+02,.2076657D-01,.3790000D+02,.0000000D+00 .,.8524012D+03,.8733411D+03,.0000000D+00,.6953030D+03,.0000000D+00 .,.9275149D+03,.1036299D+04,.0000000D+00,.7734198D+03,.0000000D+00 .,.0000000D+00,.0000000D+00,-.4820407D+00,.5269871D-01,.000000D+00 .,-.4120379D+00,.5831442D-01,.2161068D+00,.000000D+00,.5962126D+00/ c total with 1s0 np : 6177.2837686797 4761 f^{2}_{pi}= 0.074500 c total without : 5990.9917387491 4761 f^{2}_{pi}= 0.074500 DO 1 I=1,8 DO 1 J=1,5 PAR(I,J)=PARTR(J,I) 1 CONTINUE WRITE(6,*) ' Nijmegen extended soft-core potential parameters' DO 4 I=1,8 4 WRITE(6,5) (PAR(I,J),J=1,5) 5 FORMAT(5(1X,D15.7)) RETURN END ************************************************************************ SUBROUTINE NYMCPFF IMPLICIT REAL*8 (A-H,O-Z) COMMON/COPLNG/ALPV,THPV,PV1, FPI,FETA,FETP,FPI2,FETA2,FETP2, . ALVE,THV ,GV1, GRO,GOM,GFI, GRO2,GOM2,GFI2, . ALVM, FV1, FRO,FOM,FFI, FRO2,FOM2,FFI2, . ALGS,THGS,GS1, GA0,GEP,GF0, GA02,GEP2,GF02, . GFPRO,GFPOM,GFPFI, GFMRO,GFMOM,GFMFI, . FPIC2,GROC2,FROC2,GFPROC,GFMROC,GA0C2 COMMON/KAPPAV/AKRO,AKOM,AKFI COMMON/SCALIN/AMT,AMPV,AMT24 DATA SR3/1.7320508075688772935D0/ data iprint/0/ ** pseudovector coupling constants PV8 = FPI * (4D0*ALPV-1D0)/SR3 COST = DCOS(THPV) SINT = DSIN(THPV) FETA = COST*PV8 - SINT*PV1 FETP = SINT*PV8 + COST*PV1 FPI2 = FPI*FPI/(AMPV*AMPV) FETA2= FETA*FETA/(AMPV*AMPV) FETP2= FETP*FETP/(AMPV*AMPV) ** vector coupling constants GV8 = GRO * (4D0*ALVE-1D0)/SR3 COST = DCOS(THV) SINT = DSIN(THV) GFI = COST*GV8 - SINT*GV1 GOM = SINT*GV8 + COST*GV1 GRO2 = GRO*GRO GOM2 = GOM*GOM GFI2 = GFI*GFI ** tensor coupling constants FRO = AKRO*GRO FV8 = (GRO+FRO) * (4D0*ALVM-1D0)/SR3 - GV8 COST = DCOS(THV) SINT = DSIN(THV) FFI = COST*FV8 - SINT*FV1 FOM = SINT*FV8 + COST*FV1 AKOM = FOM/GOM AKFI = FFI/GFI FRO2 = FRO*FRO/AMT24 FOM2 = FOM*FOM/AMT24 FFI2 = FFI*FFI/AMT24 GFPRO = (GRO*FRO+FRO*GRO)/(2D0*AMT) GFPOM = (GOM*FOM+FOM*GOM)/(2D0*AMT) GFPFI = (GFI*FFI+FFI*GFI)/(2D0*AMT) GFMRO = (GRO*FRO-FRO*GRO)/(2D0*AMT) GFMOM = (GOM*FOM-FOM*GOM)/(2D0*AMT) GFMFI = (GFI*FFI-FFI*GFI)/(2D0*AMT) ** scalar coupling constants COST = DCOS(THGS) SINT = DSIN(THGS) GS8 = COST*GF0 + SINT*GEP GS1 =-SINT*GF0 + COST*GEP ALGS = (SR3*GS8+GA0)/(4D0*GA0) GA02 = GA0*GA0 GEP2 = GEP*GEP GF02 = GF0*GF0 if(iprint.eq.0) then iprint=1 write(*,'(a,5f13.6)') ' ga0,gep,gf0:',ga0,gep,gf0,gs8,gs1 write(*,'(a,3f13.6)') ' fpi,fet,fx0:',fpi,feta,fetp write(*,'(a,3f13.6)') ' gro,gom,gfi:',gro,gom,gfi write(*,'(a,3f13.6)') ' fro,fom,ffi:',fro,fom,ffi write(*,'(a,3f13.6)') ' kro,kom,kfi:',akro,akom,akfi endif RETURN END ************************************************************************ SUBROUTINE CHIRAL C* Imposes theoretical constraints on the pair coupling constants C* based on chiral Langrangians (see notes V. Stoks), C* using the linear sigma model. IMPLICIT REAL*8 (A-H,O-Z) COMMON/MESONM/AMPI,AMETA,AMETP, AMRO,AMOM,AMFI, AMA0,AMEP,AMF0, . AMPIC,AMROC,AMA0C,AWPIC,AWROC,AWA0C, AVSC COMMON/MEANM/ PIM,ROM,A0M COMMON/SCALIN/AMT,AMPV,AMT24 COMMON/COPLNG/ALPV,THPV,PV1, FPI,FETA,FETP,FPI2,FETA2,FETP2, . ALVE,THV ,GV1, GRO,GOM,GFI, GRO2,GOM2,GFI2, . ALVM, FV1, FRO,FOM,FFI, FRO2,FOM2,FFI2, . ALGS,THGS,GS1, GA0,GEP,GF0, GA02,GEP2,GF02, . GFPRO,GFPOM,GFPFI, GFMRO,GFMOM,GFMFI, . FPIC2,GROC2,FROC2,GFPROC,GFMROC,GA0C2 COMMON/GPAIR/ GPIETA,GPIETP,GSISI, GPISIG,GPIPOM, . GPIPI0,GPIPI1,FPIPI1,GPIRO0,GPIRO1,GPIOM COMMON/PARAMS/PAR(8,5) DATA FPIDEC/92.4D0/, GAXI/1.2573D0/ DATA PI4/12.5663706143592D0/, SQ4PI/3.54490770181103D0/ DATA ICAL/0/, SQ2/1.4142135623731D0/, SQ6/2.4494897427832D0/ DATA IC/-1/, iw/0/ SAVE AMPV2,ROM2,AMA1,A1M2,ZM1,CHIRPAR,GROCHIR IF(ICAL.EQ.0) THEN ICAL=1 AMPV2=AMPV*AMPV ROM2=ROM*ROM AMA1=1230D0 A1M2=AMA1*AMA1 ZM1 =ROM2/A1M2 CHIRPAR = ZM1*(GAXI-IC*ZM1)/(2D0*FPIDEC**2) GROCHIR = DSQRT(1D0-ZM1)*ROM/(2D0*FPIDEC) ENDIF c if(iw.eq.0) write(*,*) ' axial:',gaxi/(2d0*fpidec*sq4pi),fpi/ampv c if(iw.eq.0) write(*,*) 'vector:',grochir/sq4pi,gro GPIPI1 =-IC*AMPV2*CHIRPAR/PI4 FPIPI1 = GPIPI1 * FRO/GRO GPISIG = (A1M2/ROM2-2D0)*AMA1/ROM*AMPV2*CHIRPAR/PI4 GPIRO1 = GRO*AMPV/FPIDEC*(GAXI-IC)/SQ4PI GPIPI0 =-(AMT/FPIDEC/PI4-GEP/SQ4PI*AMA1/ROM)*AMPV/(2D0*FPIDEC) GSISI = 3D0*ZM1*GPIPI0 GPIETA = 2D0*ROM/AMA1*GPIPI0 gpieta = 0d0 if(iw.eq.0) then iw=1 write(*,'(a,2f13.6)') ' gpipi,fpipi:',gpipi1,fpipi1 write(*,'(a,2f13.6)') ' gpisi,gpiro:',gpisig,gpiro1 write(*,'(a,3f13.6)') ' p02,si2,pet:',gpipi0,gsisi,gpieta endif C* Store the pair couplings back into the parameter matrix PAR(7,1) = GPIETA PAR(7,2) = GPIETP PAR(7,3) = GSISI PAR(7,4) = GPISIG PAR(7,5) = GPIPOM PAR(8,1) = GPIPI0 PAR(8,2) = GPIPI1 PAR(8,3) = FPIPI1 PAR(8,4) = GPIRO0 PAR(8,5) = GPIRO1 PAR(6,5) = GPIOM RETURN END ************************************************************************ SUBROUTINE AMFACNN IMPLICIT REAL*8 (A-H,O-Z) COMMON/AMCOEF/ALF,REDM, AMY,AMY2,AMYI,AMYI2, AMN,AMN2,AMNI,AMNI2, . AMYPN,AMYMN, AMYPN2,AMYMN2, AY2PN2,AY2MN2, . AMYN,AMYNI,AMYNI2, AYPNI2,AYMNI2 AMY2 = AMY*AMY AMYI = 1D0/AMY AMYI2 = AMYI*AMYI AMN2 = AMN*AMN AMNI = 1D0/AMN AMNI2 = AMNI*AMNI AMYPN = AMY+AMN AMYMN = AMY-AMN AMYPN2= AMYPN*AMYPN AMYMN2= AMYMN*AMYMN AY2PN2= AMY2+AMN2 AY2MN2= AMY2-AMN2 AMYN = AMY*AMN AMYNI = AMYI*AMNI AMYNI2= AMYNI*AMYNI AYPNI2= AMYI2+AMNI2 AYMNI2= AMYI2-AMNI2 REDM = AMYN/AMYPN ALF = 4D0*REDM/AMYPN IF(AMY.EQ.AMN) ALF = 1D0 RETURN END ************************************************************************ SUBROUTINE RSPACFF(RH,TYPE,NADIA) ** pseudoscalar functions P1-4: pi+, pi0, eta, eta' ** scalar functions S1-4: a0+, a00, epsilon, f0 ** vector functions VV1-4: g*g in rho+, rho0, omega, phi ** VT1-4: g*f in rho+, rho0, omega, phi ** TT1-4: f*f in rho+, rho0, omega, phi ** VV5: g*g scalar part of rho+ IMPLICIT REAL*8 (A-H,O-Z) CHARACTER TYPE*2 LOGICAL BRDRO,BRDROC,BRDEP COMMON/MESONM/AMPI,AMETA,AMETP, AMRO,AMOM,AMFI, AMA0,AMEP,AMF0, . AMPIC,AMROC,AMA0C,AWPIC,AWROC,AWA0C, AVSC COMMON/MEANM/ PIM,ROM,A0M COMMON/YUKEFF/ARO,AMR1,AROC,AMRC1,AWRC1,BRO,AMR2,BROC,AMRC2,AWRC2, . AVSC1,AVSC2, AEPS,AME1,BEPS,AME2 COMMON/BRDRHO/GAMRO,THRRO,NRO,BRDRO,GAMROC,THRROC,NROC,BRDROC COMMON/BRDEPS/GAMEP,THREP0,THREPC,NEP,BRDEP COMMON/CUTOFF/ALAM1,ALAMP1,ALAMV1,ALAMS1, . ALAM0,ALAMP0,ALAMV0,ALAMS0 COMMON/CUTTEN/ALAMT1,ALAVT1,ALAMT0,ALAVT0 COMMON/SPCTRL/RHB,AMTHR2,AMMES,GAMMES,ALABDA,N COMMON/RVMAT/ P1(6),P2(6),P3(6),P4(6),S1(6),S2(6),S3(6),S4(6), . VV1(6),VV2(6),VV3(6),VV4(6),VT1(6),VT2(6), . VT3(6),VT4(6),TT1(6),TT2(6),TT3(6),TT4(6),VV5(6) COMMON/RBMAT/ RES(6,2),REV(6,6) COMMON/RPINA/ F3PIS(3,2),F3PIT(3,2) RHB = RH C** Pseudoscalar mesons pi, eta, eta' CALL REGGEFF(RH,AMPI,ALAMP1,P2) IF(NADIA.EQ.1) CALL RNADIA(RH,TYPE,P2,F3PIS(1,2),F3PIT(1,2)) CALL REGGEFF(RH,AMETA,ALAMP0,P3) CALL REGGEFF(RH,AMETP,ALAMP0,P4) c** Scalar mesons a0(980), epsilon, f0(975) CALL REGGEFF(RH,AMA0,ALAMS1,S2) IF(BRDEP) THEN AMMES = AMEP GAMMES= GAMEP ALABDA= ALAMS0 N = NEP AMTHR2=THREP0*THREP0 CALL BROAD(THREP0,RES(1,1),6) AMTHR2=THREPC*THREPC CALL BROAD(THREPC,RES(1,2),6) DO 10 I=1,6 S3(I) = (RES(I,1) + 2D0*RES(I,2))/3D0 10 CONTINUE ELSE CALL REGGEFF(RH,AME1,ALAMS0,RES(1,1)) CALL REGGEFF(RH,AME2,ALAMS0,RES(1,2)) DO 11 I=1,6 S3(I) = AEPS*RES(I,1) + BEPS*RES(I,2) 11 CONTINUE ENDIF CALL REGGEFF(RH,AMF0,ALAMS0,S4) C** Vector mesons rho, omega, phi IF(BRDRO) THEN AMTHR2= THRRO*THRRO AMMES = AMRO GAMMES= GAMRO N = NRO ALABDA= ALAMV1 CALL BROAD(THRRO,VV2,6) IF(ALAMT1.NE.ALAMV1) THEN ALABDA= ALAVT1 CALL BROAD(THRRO,VT2,6) ALABDA= ALAMT1 CALL BROAD(THRRO,TT2,6) ELSE DO 21 I=1,6 VT2(I) = VV2(I) TT2(I) = VV2(I) 21 CONTINUE ENDIF ELSE CALL REGGEFF(RH,AMR1,ALAMV1,REV(1,1)) CALL REGGEFF(RH,AMR2,ALAMV1,REV(1,2)) IF(ALAMT1.NE.ALAMV1) THEN CALL REGGEFF(RH,AMR1,ALAVT1,REV(1,3)) CALL REGGEFF(RH,AMR2,ALAVT1,REV(1,4)) CALL REGGEFF(RH,AMR1,ALAMT1,REV(1,5)) CALL REGGEFF(RH,AMR2,ALAMT1,REV(1,6)) ENDIF DO 22 I=1,6 VV2(I) = ARO*REV(I,1) + BRO*REV(I,2) IF(ALAMT1.NE.ALAMV1) THEN VT2(I) = ARO*REV(I,3) + BRO*REV(I,4) TT2(I) = ARO*REV(I,5) + BRO*REV(I,6) ELSE VT2(I) = ARO*REV(I,1) + BRO*REV(I,2) TT2(I) = ARO*REV(I,1) + BRO*REV(I,2) ENDIF 22 CONTINUE ENDIF CALL REGGEFF(RH,AMOM,ALAMV0,VV3) CALL REGGEFF(RH,AMFI,ALAMV0,VV4) IF(ALAMT0.NE.ALAMV0) THEN CALL REGGEFF(RH,AMOM,ALAVT0,VT3) CALL REGGEFF(RH,AMOM,ALAMT0,TT3) CALL REGGEFF(RH,AMFI,ALAVT0,VT4) CALL REGGEFF(RH,AMFI,ALAMT0,TT4) ELSE DO 23 I=1,6 VT3(I) = VV3(I) TT3(I) = VV3(I) VT4(I) = VV4(I) TT4(I) = VV4(I) 23 CONTINUE ENDIF C** Isospin-1 mesons (pi+, rho+, a0+) IF(TYPE.EQ.'NP' .OR. TYPE.EQ.'PN') THEN CALL REGGEFF(RH,AWPIC,ALAMP1,P1) IF(NADIA.EQ.1) CALL RNADIA(RH,TYPE,P1,F3PIS(1,1),F3PIT(1,1)) IF(BRDROC) THEN AMTHR2= THRROC*THRROC AMMES = AMROC GAMMES= GAMROC N = NROC ALABDA= ALAMV1 CALL BROAD(THRROC,VV1,6) IF(ALAMT1.NE.ALAMV1) THEN ALABDA= ALAVT1 CALL BROAD(THRROC,VT1,6) ALABDA= ALAMT1 CALL BROAD(THRROC,TT1,6) ELSE DO 31 I=1,6 VT1(I) = VV1(I) TT1(I) = VV1(I) 31 CONTINUE ENDIF DO 32 I=1,6 VV5(I) = AVSC*VV1(I) 32 CONTINUE ELSE DO 33 I=1,6 VV1(I) = AROC*REV(I,1) + BROC*REV(I,2) IF(ALAMT1.NE.ALAMV1) THEN VT1(I) = AROC*REV(I,3) + BROC*REV(I,4) TT1(I) = AROC*REV(I,5) + BROC*REV(I,6) ELSE VT1(I) = AROC*REV(I,1) + BROC*REV(I,2) TT1(I) = AROC*REV(I,1) + BROC*REV(I,2) ENDIF VV5(I) = AVSC1*AROC*REV(I,1) + AVSC2*BROC*REV(I,2) 33 CONTINUE ENDIF CALL REGGEFF(RH,AWA0C,ALAMS1,S1) ENDIF RETURN END ************************************************************************ SUBROUTINE REGGEFF(RH,AM,ALAM,FIFUNS) IMPLICIT REAL*8 (A-H,O-Z) DIMENSION FIFUNS(6) DATA ZMATCH/1D-2/, SQPI/1.7724538509055160273D0/ RATM=AM/ALAM RATM2=RATM*RATM RATM3=RATM*RATM2 EXPM=FDEXP(RATM2) AM3=AM*AM*AM AM5=AM*AM*AM3 Z=RH*ALAM IF(Z.GT.ZMATCH) THEN XLAM=0.5D0*RH*ALAM XLAM2=XLAM*XLAM ELAM=FDEXP(-XLAM2) X=AM*RH ERFCM=FDERFC(RATM-XLAM) ERFCP=FDERFC(RATM+XLAM) EPX=0D0 IF(ERFCP.NE.0) EPX=FDEXP(X) EMX=FDEXP(-X) PHI0C = EXPM*(EMX*ERFCM-EPX*ERFCP)/(2D0*X) PHI0L =(EXPM*(EMX*ERFCM*(1D0+X)-EPX*ERFCP*(1D0-X))- . 4D0/SQPI*XLAM*ELAM ) / (2D0*X*X*X) PHI0T =(EXPM*(EMX*ERFCM*(1D0+X+X*X/3D0) . -EPX*ERFCP*(1D0-X+X*X/3D0)) . -8D0/SQPI*XLAM*(0.5D0+XLAM2/3D0)*ELAM) / (2D0*X*X*X) DEL1C =-ELAM/(RATM3*2D0*SQPI) DEL1L =-ELAM/(RATM2*RATM3*4D0*SQPI) DEL1T =-ELAM*XLAM2/(RATM2*RATM3*6D0*SQPI) DEL2C = 3D0*(DEL1T-DEL1L) ELSE A=RATM B=-FDERFC(A)*A*EXPM*SQPI A2=A*A A3=A2*A A4=A2*A2 A5=A3*A2 A6=A4*A2 A8=A6*A2 A10=A8*A2 Z2=Z*Z Z4=Z2*Z2 Z6=Z4*Z2 Z8=Z6*Z2 PHI0C = (16*A8*B+16*A8-8*A6+12*A4-30*A2+105)*Z8/5806080D0+ . (8*A6*B+8*A6-4*A4+6*A2-15)*Z6/40320D0+ . (4*A4*B+4*A4-2*A2+3)*Z4/480D0+ . (2*A2*B+2*A2-1)*Z2/12D0+B+1D0 PHI0C = PHI0C/(SQPI*A) PHI0T = (32*A10*B+32*A10-16*A8+24*A6-60*A4+210*A2-945) . *Z8/47900160D0+ . (16*A8*B+16*A8-8*A6+12*A4-30*A2+105)*Z6/362880D0+ . (8*A6*B+8*A6-4*A4+6*A2-15)*Z4/5040D0+ . (4*A4*B+4*A4-2*A2+3)*Z2/180D0 PHI0T = PHI0T/(SQPI*A3) PHI0L =-(32*A10*B+32*A10-16*A8+24*A6-60*A4+210*A2-945) . *Z8/127733760D0- . (16*A8*B+16*A8-8*A6+12*A4-30*A2+105)*Z6/725760D0- . (8*A6*B+8*A6-4*A4+6*A2-15)*Z4/6720D0- . (4*A4*B+4*A4-2*A2+3)*Z2/120D0-(2*A2*B+2*A2-1D0)/6D0 PHI0L = PHI0L/(SQPI*A3) DEL1C =-Z8/12288D0+Z6/768D0-Z4/64D0+Z2/8D0-0.5D0 DEL1C = DEL1C/(SQPI*A3) DEL1T = Z8/9216D0-Z6/768D0+Z4/96D0-Z2/24D0 DEL1T = DEL1T/(SQPI*A5) DEL1L =-Z8/24576D0+Z6/1536D0-Z4/128D0+Z2/16D0-0.25D0 DEL1L = DEL1L/(SQPI*A5) DEL2C = 3D0*(DEL1T-DEL1L) ENDIF FIFUNS(1) = AM*PHI0C FIFUNS(2) = AM3*(PHI0C+DEL1C) FIFUNS(3) = AM5*(PHI0C+DEL1C+DEL2C) FIFUNS(4) = AM3*PHI0L FIFUNS(5) = AM5*(PHI0L+DEL1L) FIFUNS(6) = AM3*PHI0T c FIFUNS(7) = AM5*(PHI0T+DEL1T) RETURN END ************************************************************************ SUBROUTINE BROAD(AMTHR,SUM,NDIM) ** Calculates integrals for broad mesons using Gauss integration ** on (mthr,infinity) by mapping (mthr,3*mthr) on (-1,1) with 16 ** points, and (3*mthr,20*mthr) on (-1,1) with 8 points. IMPLICIT REAL*8(A-H,O-Z) DIMENSION SUM(NDIM),FUNCW(6),YB(3),NP(2),X(20,2),W(20,2) DATA ICAL/0/, YB/1D0,3D0,20D0/, NP/16,8/, NI/2/ SAVE X,W IF(ICAL.EQ.0) THEN CALL GAUSS(NP(1),X(1,1),W(1,1),1) CALL GAUSS(NP(2),X(1,2),W(1,2),1) ICAL=1 ENDIF DO 1 I=1,NDIM SUM(I)=0D0 1 CONTINUE DO 2 IN=1,NI BMA=0.5D0*(YB(IN+1)-YB(IN))*AMTHR BPA=0.5D0*(YB(IN+1)+YB(IN))*AMTHR DO 3 IY=1,NP(IN) AMY=BMA*X(IY,IN)+BPA CALL BRDFUN(FUNCW,AMY) DO 4 K=1,NDIM SUM(K)=SUM(K)+FUNCW(K)*W(IY,IN)*BMA 4 CONTINUE 3 CONTINUE 2 CONTINUE RETURN END ************************************************************************ SUBROUTINE BRDFUN(FINTGR,AMP) IMPLICIT REAL*8 (A-H,O-Z) DIMENSION FI(6),FINTGR(6) COMMON/SPCTRL/RB,AMTHR2,AMMES,GAMMES,ALAM,N DATA PI/3.14159265358979323846D0/ AMP2=AMP*AMP AMMES2=AMMES*AMMES FACG=AMMES*GAMMES/(AMMES2-AMTHR2)**(N+0.5D0) SPEC=FACG*(AMP2-AMTHR2)**(N+0.5D0)/ . ( (AMP2-AMMES2)**2+FACG*FACG*(AMP2/AMMES2)**(2D0*N)* . (AMP2-AMTHR2)**(2D0*N+1D0) ) / PI CALL REGGEFF(RB,AMP,ALAM,FI) DO 1 I=1,6 FINTGR(I)=2D0*AMP*SPEC*FI(I) 1 CONTINUE RETURN END ************************************************************************ SUBROUTINE RNADIA(RH,TYPE,F2PI,F3PIS,F3PIT) C* Calculates the I_{3} functions with first and second derivatives C* of spin-spin and tensor OPE functions. C* These derivatives can be written in terms of C* PHI^{1}_{C}, PHI^{2}_{C}, PHI^{1}_{SO}, PHI^{0}_{T} C----------------------------------------------------------------------- IMPLICIT REAL*8 (A-H,O-Z) REAL*8 F2PI(6),FILAM(6), X(20),W(20), SUM(4), F3PIS(3),F3PIT(3) CHARACTER TYPE*2 COMMON/MESONM/AMPI,AMETA,AMETP, AMRO,AMOM,AMFI, AMA0,AMEP,AMF0, . AMPIC,AMROC,AMA0C,AWPIC,AWROC,AWA0C, AVSC COMMON/MEANM/ PIM,ROM,A0M COMMON/CUTOFF/ALAM1,ALAMP1,ALAMV1,ALAMS1, . ALAM0,ALAMP0,ALAMV0,ALAMS0 DATA NP/12/, ICAL/0/, TWONPI/0.636619772367581343D0/ IF(ICAL.EQ.0) THEN CALL GAUSS(NP,X,W,1) DO 1 I=1,NP W(I)=2D0*W(I)/(1D0-X(I))**2 X(I)=(1D0+X(I))/(1D0-X(I)) 1 CONTINUE ICAL=1 ENDIF IF(TYPE.EQ.'PP' .OR. TYPE.EQ.'NN') THEN AMES=AMPI ELSEIF(TYPE.EQ.'NP' .OR. TYPE.EQ.'PN') THEN AMES=AWPIC ENDIF DO 10 K=1,4 10 SUM(K)=0D0 DO 11 IY=1,NP PLAM=X(IY)*PIM PLAM2=PLAM*PLAM CALL PHIFUN(RH,AMES,ALAMP1,PLAM,FILAM) SUM(1)=SUM(1)+(F2PI(2)-FILAM(2))/PLAM2*W(IY) SUM(2)=SUM(2)+(F2PI(3)-FILAM(3))/PLAM2*W(IY) SUM(3)=SUM(3)+(F2PI(5)-FILAM(5))/PLAM2*W(IY) SUM(4)=SUM(4)+(F2PI(6)-FILAM(6))/PLAM2*W(IY) 11 CONTINUE ** spin (0th, 1st, 2nd) and tensor (0th, 1st, 2nd) F3PIS(1) = PIM*TWONPI * SUM(1) F3PIS(2) =-PIM*TWONPI * SUM(3)*RH F3PIS(3) = PIM*TWONPI *(SUM(2)+2D0*SUM(3)) F3PIT(1) = PIM*TWONPI * SUM(4) F3PIT(2) =-PIM*TWONPI *(3D0/RH*SUM(4)+RH/3D0*SUM(3)) F3PIT(3) = PIM*TWONPI *(12D0/(RH*RH)*SUM(4) . +(5D0*SUM(3)+SUM(2))/3D0) RETURN END *********************************************************************** SUBROUTINE PHIFUN(RH,AM,ALAM,PLAM,PHI) ** Basic function phi^{n}_{X} stored as ** PHI = (m phi^{0}_{C}, m^{3} phi^{1}_{C}, m^{5} phi^{2}_{C}, ** m^{3} phi^{0}_{SO}, m^{5} phi^{1}_{SO}, ** m^{3} phi^{0}_{T}, m^{5} phi^{1}_{T}) IMPLICIT REAL*8 (A-H,O-Z) DIMENSION PHI(6) DATA SRPI/1.7724538509055160273D0/ AMP = DSQRT(AM*AM+PLAM*PLAM) AMP2 = AMP*AMP AMP3 = AMP*AMP2 RATM = AMP/ALAM RATM2= RATM*RATM RATM3= RATM*RATM2 EXPM = FDEXP(AM*AM/(ALAM*ALAM)) EXPL = FDEXP(-PLAM*PLAM/(ALAM*ALAM)) XLAM = 0.5D0*RH*ALAM XLAM2= XLAM*XLAM ELAM = FDEXP(-XLAM2) X = AMP*RH X2= X*X X3= X*X2 ERFCM = FDERFC(RATM-XLAM) ERFCP = FDERFC(RATM+XLAM) EPX = 0D0 IF(ERFCP.NE.0) EPX = FDEXP(X) EMX = FDEXP(-X) PHI0C = EXPM*(EMX*ERFCM-EPX*ERFCP)/(2D0*X) PHI0L = (EXPM*(EMX*ERFCM*(1D0+X)-EPX*ERFCP*(1D0-X)) . - 4D0/SRPI*XLAM*ELAM*EXPL ) / (2D0*X3) PHI0T =(EXPM*(EMX*ERFCM*(1D0+X+X2/3D0)-EPX*ERFCP*(1D0-X+X2/3D0)) . - 8D0/SRPI*XLAM*(0.5D0+XLAM2/3D0)*ELAM*EXPL ) / (2D0*X3) DEL1C =-ELAM/(RATM3*2D0*SRPI) * EXPL PHI1C = PHI0C+DEL1C DEL1L =-ELAM/(RATM2*RATM3*4D0*SRPI) * EXPL PHI1L = PHI0L+DEL1L DEL1T =-ELAM*XLAM2/(RATM2*RATM3*6D0*SRPI) * EXPL PHI1T = PHI0T+DEL1T PHI2C = PHI1C+3D0*(DEL1T-DEL1L) PHI(1) = AMP *PHI0C PHI(2) = AMP3*PHI1C PHI(3) = AMP3*PHI2C*AMP2 PHI(4) = AMP3*PHI0L PHI(5) = AMP3*PHI1L*AMP2 PHI(6) = AMP3*PHI0T c PHI(7) = AMP3*PHI1T*AMP2 RETURN END ************************************************************************ SUBROUTINE POTFNN(RH,TYPE,ISIGN,FACISO,NLOC,NADIA) IMPLICIT REAL*8 (A-H,O-Z) CHARACTER TYPE*2 DIMENSION DUM1(9),DUM2(9) COMMON/POT2/ VC,VSS,VT,VLS,VLSA,VQ12,FIC0,FIC1,FIC2 COMMON/VOPE/ VPIS,VPIT,VPIQ,FPIS0,FPIS1,FPIS2,FPIT0,FPIT1,FPIT2 COMMON/VENER/ PLMEV2,VEPIS,VEPIT COMMON/COPLNG/ALPV,THPV,PV1, FPI,FETA,FETP,FPI2,FETA2,FETP2, . ALVE,THV ,GV1, GRO,GOM,GFI, GRO2,GOM2,GFI2, . ALVM, FV1, FRO,FOM,FFI, FRO2,FOM2,FFI2, . ALGS,THGS,GS1, GA0,GEP,GF0, GA02,GEP2,GF02, . GFPRO,GFPOM,GFPFI, GFMRO,GFMOM,GFMFI, . FPIC2,GROC2,FROC2,GFPROC,GFMROC,GA0C2 COMMON/SCALIN/AMT,AMPV,AMT24 COMMON/POMRON/GPOM,FPOM2,AMPOM,AMPOM2,AMPOM4, . GA2D,FA2D2,AMA2D,AMA2D2,AMA2D4 COMMON/AMCOEF/ALF,REDM, AMY,AMY2,AMYI,AMYI2, AMN,AMN2,AMNI,AMNI2, . AMYPN,AMYMN, AMYPN2,AMYMN2, AY2PN2,AY2MN2, . AMYN,AMYNI,AMYNI2, AYPNI2,AYMNI2 COMMON/RVMAT/ P1(6),P2(6),P3(6),P4(6),S1(6),S2(6),S3(6),S4(6), . VV1(6),VV2(6),VV3(6),VV4(6),VT1(6),VT2(6), . VT3(6),VT4(6),TT1(6),TT2(6),TT3(6),TT4(6),VV5(6) COMMON/RPINA/ F3PIS(3,2),F3PIT(3,2) EQUIVALENCE (VC,DUM1(1)), (VPIS,DUM2(1)) DO 1 ID=1,9 DUM1(ID)=0D0 DUM2(ID)=0D0 1 CONTINUE RH2=RH*RH C* Pseudoscalar-meson potential F2PION=FPI2 IF(NADIA.EQ.1) THEN ** Nonadiabatic correction to OPE, energy multiplication outside FISP0 =-0.5D0*FPI2*F3PIS(1,2)/3D0 FISP1 =-0.5D0*FPI2*F3PIS(2,2)/3D0 FISP2 =-0.5D0*FPI2*F3PIS(3,2)/3D0 FITP0 =-0.5D0*FPI2*F3PIT(1,2) FITP1 =-0.5D0*FPI2*F3PIT(2,2) FITP2 =-0.5D0*FPI2*F3PIT(3,2) VEPIS = 0.5D0*FPI2*F3PIS(1,2)/(3D0*REDM) VEPIT = 0.5D0*FPI2*F3PIT(1,2)/REDM ELSEIF(NLOC.EQ.1) THEN ** Energy-dependent OPE, added separately CONTINUE ELSEIF(NLOC.EQ.2) THEN ** Nonlocal spin-spin and tensor FFNL = REDM*FPI2*AMYNI/2D0 FPIS0=-FFNL*P2(2)/3D0 FPIS1= FFNL*P2(5)*RH/3D0 FPIS2=-FFNL*(P2(3)+2D0*P2(5))/3D0 FPIT0=-FFNL*P2(6) FPIT1= FFNL*(3D0*P2(6)/RH+P2(5)*RH/3D0) FPIT2=-FFNL*(12D0*P2(6)/RH2+(5D0*P2(5)+P2(3))/3D0) ENDIF VPIS = FPI2*P2(2)/3D0 VPIT = FPI2*P2(6) VPVS = (FETA2*P3(2)+FETP2*P4(2))/3D0 VPVT = (FETA2*P3(6)+FETP2*P4(6)) FPI2=F2PION C* Scalar-meson potential VSC =-(GA02*S2(1)+GEP2*S3(1)+GF02*S4(1)) + . (GA02*S2(2)+GEP2*S3(2)+GF02*S4(2)) * AYPNI2/8D0 VSLS =-(GA02*S2(4)+GEP2*S3(4)+GF02*S4(4)) * AYPNI2/4D0 VSLA =-(GA02*S2(4)+GEP2*S3(4)+GF02*S4(4)) * AYMNI2/4D0 VSQ12=-(GA02*S2(6)+GEP2*S3(6)+GF02*S4(6)) * AMYNI2/16D0 . *3D0/RH2 FICS0= REDM*(GA02*S2(1)+GEP2*S3(1)+GF02*S4(1))*AMYNI/2D0 FICS1=-REDM*(GA02*S2(4)+GEP2*S3(4)+GF02*S4(4))*AMYNI/2D0*RH FICS2= REDM*(GA02*S2(2)+GEP2*S3(2)+GF02*S4(2))*AMYNI/2D0 . + REDM*(GA02*S2(4)+GEP2*S3(4)+GF02*S4(4))*AMYNI C* Diffractive contribution XI0 = AMPOM*RH XI02 = XI0*XI0 EXI0 = 0D0 IF(XI02.LT.170.) EXI0 = FDEXP(-XI02) XI1 = AMA2D*RH XI12 = XI1*XI1 EXI1 = 0D0 IF(XI12.LT.170.) EXI1 = FDEXP(-XI12) VDC = FPOM2 * (1D0+AMPOM2*(1.5D0-XI02)*AYPNI2/2D0) * EXI0 + . FA2D2 * (1D0+AMA2D2*(1.5D0-XI12)*AYPNI2/2D0) * EXI1*ISIGN VDLS = (FPOM2 * AMPOM2 * EXI0 + . FA2D2 * AMA2D2 * EXI1*ISIGN) * AYPNI2/2D0 VDLA = (FPOM2 * AMPOM2 * EXI0 + . FA2D2 * AMA2D2 * EXI1*ISIGN) * AYMNI2/2D0 VDQ12= (FPOM2 * AMPOM4 * EXI0 + . FA2D2 * AMA2D4 * EXI1*ISIGN) * AMYNI2/4D0 FICD0=-REDM*(FPOM2*EXI0 + ISIGN*FA2D2*EXI1) * AMYNI/2D0 FICD1= REDM*(FPOM2*EXI0*XI02 + ISIGN*FA2D2*EXI1*XI12) * AMYNI/RH FICD2= REDM*(FPOM2*EXI0*(1D0-2D0*XI02)*AMPOM2 + . ISIGN*FA2D2*EXI1*(1D0-2D0*XI12)*AMA2D2) * AMYNI C* Vector-meson potential VVC = (GRO2*VV2(1)+GOM2*VV3(1)+GFI2*VV4(1)) . +(GRO2*VV2(2)+GOM2*VV3(2)+GFI2*VV4(2)) * AMYNI/(2D0*ALF) . +(gro2*vv2(3)+gom2*vv3(3)+gfi2*vv4(3)) * amyni2/64d0 . +(GFPRO*VT2(2)+GFPOM*VT3(2)+GFPFI*VT4(2)) / (4D0*REDM) . +(gfpro*vt2(3)+gfpom*vt3(3)+gfpfi*vt4(3))*amyni/(16d0*redm) . +(GFMRO*VT2(2)+GFMOM*VT3(2)+GFMFI*VT4(2)) * AMYMN*AMYNI/4D0 . +(FRO2*TT2(3)+FOM2*TT3(3)+FFI2*TT4(3))*AMYNI/4D0 VVS = (GRO2*VV2(2)+GOM2*VV3(2)+GFI2*VV4(2)) * AMYNI/6D0 . +(GFPRO*VT2(2)+GFPOM*VT3(2)+GFPFI*VT4(2))*AMYPN*AMYNI/6D0 . -(GFMRO*VT2(2)+GFMOM*VT3(2)+GFMFI*VT4(2))*AMYMN*AMYNI/6D0 . +(FRO2*TT2(2)+FOM2*TT3(2)+FFI2*TT4(2))*2D0/3D0 C . +(FRO2*TT2(3)+FOM2*TT2(3)+FFI2*TT2(3))*AMYNI/12D0 VVT =-(GRO2*VV2(6)+GOM2*VV3(6)+GFI2*VV4(6)) * AMYNI/4D0 . -(GFPRO*VT2(6)+GFPOM*VT3(6)+GFPFI*VT4(6))*AMYPN*AMYNI/4D0 . +(GFMRO*VT2(6)+GFMOM*VT3(6)+GFMFI*VT4(6))*AMYMN*AMYNI/4D0 . -(FRO2*TT2(6)+FOM2*TT3(6)+FFI2*TT4(6)) C . -(FRO2*TT2(7)+FOM2*TT2(7)+FFI2*TT2(7))*AMYNI/8D0 VVLS =-(GRO2*VV2(4)+GOM2*VV3(4)+GFI2*VV4(4))*(.5D0+1D0/ALF)*AMYNI . -(gro2*vv2(5)+gom2*vv3(5)+gfi2*vv4(5))*amyni2/16d0 . -(GFPRO*VT2(4)+GFPOM*VT3(4)+GFPFI*VT4(4))/REDM . -(gfpro*vt2(5)+gfpom*vt3(5)+gfpfi*vt4(5))*amyni/(4d0*redm) . -(FRO2*TT2(5)+FOM2*TT3(5)+FFI2*TT4(5))*(.5D0+1D0/ALF)*AMYNI VVLA =-(GRO2*VV2(4)+GOM2*VV3(4)+GFI2*VV4(4))*AYMNI2/4D0 . +(GFMRO*VT2(4)+GFMOM*VT3(4)+GFMFI*VT4(4))/REDM . +(FRO2*TT2(5)+FOM2*TT3(5)+FFI2*TT4(5))*AYMNI2/4D0 VVQ12=((GRO2*VV2(6)+GOM2*VV3(6)+GFI2*VV4(6)) . +(GFPRO*VT2(6)+GFPOM*VT3(6)+GFPFI*VT4(6))*4D0*AMYPN . +(FRO2*TT2(6)+FOM2*TT3(6)+FFI2*TT4(6))*8D0*AMYPN2) . *AMYNI2/16D0 * 3D0/RH2 FICV0= REDM*(GRO2*VV2(1)+GOM2*VV3(1)+GFI2*VV4(1)) * . (AYPNI2+AMYNI)/2D0 FICV1=-REDM*(GRO2*VV2(4)+GOM2*VV3(4)+GFI2*VV4(4)) * . (AYPNI2+AMYNI)/2D0 * RH FICV2= REDM*(GRO2*VV2(2)+GOM2*VV3(2)+GFI2*VV4(2)) * . (AYPNI2+AMYNI)/2D0 . + REDM*(GRO2*VV2(4)+GOM2*VV3(4)+GFI2*VV4(4)) *(AYPNI2+AMYNI) IF(TYPE.EQ.'NP' .OR. TYPE.EQ.'PN') THEN ** Pseudoscalar: pi+ F2PION=FPIC2 IF(NADIA.EQ.1) THEN ** Nonadiabatic correction to OPE, energy multiplication outside FPIS0=FPIS0 - 0.5D0*FPIC2*F3PIS(1,1)/3D0 FPIS1=FPIS1 - 0.5D0*FPIC2*F3PIS(2,1)/3D0 FPIS2=FPIS2 - 0.5D0*FPIC2*F3PIS(3,1)/3D0 FPIT0=FPIT0 - 0.5D0*FPIC2*F3PIT(1,1) FPIT1=FPIT1 - 0.5D0*FPIC2*F3PIT(2,1) FPIT2=FPIT2 - 0.5D0*FPIC2*F3PIT(3,1) VEPIS=VEPIS + 0.5D0*FPIC2*F3PIS(1,1)/(3D0*REDM) VEPIT=VEPIT + 0.5D0*FPIC2*F3PIT(1,1)/REDM ELSEIF(NLOC.EQ.1) THEN ** Energy-dependent OPE, energy multiplication outside CONTINUE ELSEIF(NLOC.EQ.2) THEN FFNL = REDM*FPIC2*(AMYNI-AYPNI2/4D0) FPIS0=FPIS0 - FFNL*P1(2)/3D0 FPIS1=FPIS1 + FFNL*P1(5)*RH/3D0 FPIS2=FPIS2 - FFNL*(P1(3)+2D0*P1(5))/3D0 FPIT0=FPIT0 - FFNL*P1(6) FPIT1=FPIT1 + FFNL*(3D0*P1(6)/RH+P1(5)*RH/3D0) FPIT2=FPIT2 - FFNL*(12D0*P1(6)/RH2+(5D0*P1(5)+P1(3))/3D0) ENDIF VPIS = VPIS + FPIC2*P1(2)/3D0 VPIT = VPIT + FPIC2*P1(6) c ENDIF FPIC2=F2PION ** Scalar: a0+ VSC = VSC - (GA0C2*S1(1) - GA0C2*S1(2)*AMYNI/4D0) VSLS = VSLS - GA0C2*S1(4)*AMYNI/2D0 VSQ12= VSQ12- GA0C2*S1(6)*AMYNI2/16D0 * 3D0/RH2 FICS0= FICS0- REDM*GA0C2*S1(1) * (AYPNI2/4D0-AMYNI) FICS1= FICS1+ REDM*GA0C2*S1(4) * (AYPNI2/4D0-AMYNI)*RH FICS2= FICS2- REDM*GA0C2*(AYPNI2/4D0-AMYNI)*(S1(2)+2D0*S1(4)) ** Diffractive: pomeron (A2) VDC = VDC + FACISO*FA2D2* . (1D0+AMA2D2*(1.5D0-XI12)*AMYNI) * EXI1 VDLS = VDLS + FACISO*FA2D2 * AMA2D2*AMYNI * EXI1 VDQ12= VDQ12+ FACISO*FA2D2 * AMA2D4*AMYNI2/4D0 * EXI1 FICD0= FICD0+ REDM * FACISO*FA2D2*(AYPNI2/4D0-AMYNI)*EXI1 FICD1= FICD1- REDM * FACISO*FA2D2*(AYPNI2/4D0-AMYNI)*EXI1 . * 2D0*XI12/RH FICD2= FICD2- REDM * FACISO*FA2D2*(AYPNI2/4D0-AMYNI)*EXI1 * . AMA2D2*(2D0-4D0*XI12) ** Vector: rho+ VVC = VVC + GROC2*(VV1(1)+AMYNI*VV1(2)/(2D0*ALF)) . + groc2*amyni2/64d0*vv1(3) . + GFPROC*AMYPN*AYPNI2/8D0*VT1(2) . + gfproc*amypn*amyni2*vt1(3)/16d0 . + FROC2*AMYMN2*AYPNI2/8D0*TT1(2) . + froc2*aypni2/8d0*tt1(3) VVS = VVS + (1D0/ALF-1D0/3D0)*AMYNI/4D0 * . (GROC2*VV1(2)+AMYPN*GFPROC*VT1(2)+AMYPN2*FROC2*TT1(2)) C . + FROC2*AYPNI2/24D0*TT1(3) VVT = VVT - AMYNI/4D0 * . (GROC2*VV1(6)+AMYPN*GFPROC*VT1(6)+AMYPN2*FROC2*TT1(6)) C . - FROC2*AYPNI2/16D0*TT1(7) VVLS = VVLS - GROC2/2D0*(AYPNI2+AMYNI)*VV1(4) . - groc2*amyni2/16d0*vv1(5) . - GFPROC/2D0*AMYPN*AYPNI2*VT1(4) . - gfproc/4d0*amypn*amyni2*vt1(5) . - FROC2/2D0*AY2MN2**2*AMYNI2*TT1(4) . - froc2/2d0*tt1(5)*(aypni2+amyni) VVLA = VVLA + GFMROC/REDM*VT1(4) VVQ12= VVQ12 + AMYNI2/16D0 * 3D0/RH2 * (GROC2*VV1(6) . + 4D0*AMYPN*GFPROC*VT1(6)+8D0*AMYPN2*FROC2*TT1(6)) FICV0= FICV0 + REDM*GROC2*(AYPNI2+AMYNI)/2D0*VV1(1) FICV1= FICV1 - REDM*GROC2*(AYPNI2+AMYNI)/2D0*VV1(4)*RH FICV2= FICV2 + REDM*GROC2*(AYPNI2+AMYNI)/2D0 . * (VV1(2)+2D0*VV1(4)) C Second part of vector-vector potential (same as scalar) COUPL= GROC2*AMYMN2 VVC = VVC - COUPL*(VV5(1) - VV5(2)*AMYNI/4D0) VVLS = VVLS - COUPL*VV5(4)*AMYNI/2D0 VVQ12= VVQ12- COUPL*VV5(6)*AMYNI2/16D0 * 3D0/RH2 FICV0= FICV0- REDM*COUPL*(AYPNI2/4D0-AMYNI)*VV5(1) FICV1= FICV1+ REDM*COUPL*(AYPNI2/4D0-AMYNI)*VV5(4)*RH FICV2= FICV2- REDM*COUPL*(AYPNI2/4D0-AMYNI) . * (VV5(2)+2D0*VV5(4)) ENDIF VC = VSC + VVC + VDC VSS = VPVS + VVS VT = VPVT + VVT VLS = VSLS + VVLS + VDLS VLSA = VSLA + VVLA + VDLA VQ12 = VSQ12 + VVQ12 + VDQ12 FIC0 = FICS0 + FICV0 + FICD0 FIC1 = FICS1 + FICV1 + FICD1 FIC2 = FICS2 + FICV2 + FICD2 RETURN END *********************************************************************** FUNCTION FDEXP(X) IMPLICIT REAL*8(A-Z) IF(X.LE.-100D0) THEN FDEXP=0D0 ELSE FDEXP=DEXP(X) ENDIF RETURN END *********************************************************************** FUNCTION FDERFC(X) IMPLICIT REAL*8(A-Z) IF(X.GE. 10D0) THEN FDERFC=0D0 ELSEIF(X.LE.-10D0) THEN FDERFC=2D0 ELSE FDERFC=DERFC(X) ENDIF RETURN END *********************************************************************** SUBROUTINE DFI2DR(ND,RH,AM,ALAM,PLAM,D0FI,DFI) ** Basic function I_{2} = m phi^{0}_{C} and its ND derivatives. ** Each function is scaled by the charged pion mass squared (PIOMS2) IMPLICIT REAL*8 (A-H,O-Z) DIMENSION DFI(*) COMMON/AMSCAL/PIOMS2 DATA ZLIM/1D-2/, RATLIM/1D0/, SQPI/1.7724538509055160273D0/ AMP =DSQRT(AM*AM+PLAM*PLAM) AMP2 =AMP*AMP AMP3 =AMP*AMP2 RATM =AMP/ALAM RATM2=RATM*RATM RATM3=RATM*RATM2 RATM5=RATM2*RATM3 EXPM =FDEXP(AM*AM/(ALAM*ALAM)) EXPL =FDEXP(-PLAM*PLAM/(ALAM*ALAM)) Z=RH*ALAM IF(Z.GT.ZLIM .OR. RATM.GT.RATLIM) THEN XLAM = 0.5D0*RH*ALAM XLAM2= XLAM*XLAM ELAM = FDEXP(-XLAM2) X = AMP*RH X2= X*X X3= X*X2 ERFCM = FDERFC(RATM-XLAM) ERFCP = FDERFC(RATM+XLAM) EPX = 0D0 IF(ERFCP.NE.0) EPX = FDEXP(X) EMX = FDEXP(-X) PHI0C=EXPM*(EMX*ERFCM-EPX*ERFCP)/(2D0*X) PHI0L=(EXPM*(EMX*ERFCM*(1D0+X)-EPX*ERFCP*(1D0-X)) . - 4D0/SQPI*XLAM*ELAM*EXPL ) / (2D0*X3) PHI0T=(EXPM*(EMX*ERFCM*(1D0+X+X2/3D0)-EPX*ERFCP*(1D0-X+X2/3D0)) . - 8D0/SQPI*XLAM*(0.5D0+XLAM2/3D0)*ELAM*EXPL ) / (2D0*X3) DEL1C=-ELAM/(RATM3*2D0*SQPI) * EXPL DEL1L=-ELAM/(RATM5*4D0*SQPI) * EXPL DEL1T=-ELAM*XLAM2/(RATM5*6D0*SQPI) * EXPL ELSE A=RATM B=-FDERFC(A)*A*EXPM*SQPI A2=A*A A3=A2*A A4=A2*A2 A5=A3*A2 A6=A4*A2 A8=A6*A2 A10=A8*A2 Z2=Z*Z Z4=Z2*Z2 Z6=Z4*Z2 Z8=Z6*Z2 PHI0C = (16*A8*B+16*A8-8*A6+12*A4-30*A2+105)*Z8/5806080D0+ . (8*A6*B+8*A6-4*A4+6*A2-15)*Z6/40320D0+ . (4*A4*B+4*A4-2*A2+3)*Z4/480D0+ . (2*A2*B+2*A2-1)*Z2/12D0+B+1D0 PHI0C = PHI0C/(SQPI*A) * EXPL PHI0T = (32*A10*B+32*A10-16*A8+24*A6-60*A4+210*A2-945) . *Z8/47900160D0+ . (16*A8*B+16*A8-8*A6+12*A4-30*A2+105)*Z6/362880D0+ . (8*A6*B+8*A6-4*A4+6*A2-15)*Z4/5040D0+ . (4*A4*B+4*A4-2*A2+3)*Z2/180D0 PHI0T = PHI0T/(SQPI*A3) * EXPL PHI0L =-(32*A10*B+32*A10-16*A8+24*A6-60*A4+210*A2-945) . *Z8/127733760D0- . (16*A8*B+16*A8-8*A6+12*A4-30*A2+105)*Z6/725760D0- . (8*A6*B+8*A6-4*A4+6*A2-15)*Z4/6720D0- . (4*A4*B+4*A4-2*A2+3)*Z2/120D0-(2*A2*B+2*A2-1D0)/6D0 PHI0L = PHI0L/(SQPI*A3) * EXPL DEL1C =-Z8/12288D0+Z6/768D0-Z4/64D0+Z2/8D0-0.5D0 DEL1C = DEL1C/(SQPI*A3) * EXPL DEL1T = Z8/9216D0-Z6/768D0+Z4/96D0-Z2/24D0 DEL1T = DEL1T/(SQPI*A5) * EXPL DEL1L =-Z8/24576D0+Z6/1536D0-Z4/128D0+Z2/16D0-0.25D0 DEL1L = DEL1L/(SQPI*A5) * EXPL ENDIF D0FI = AMP*PHI0C / PIOMS2 DFI(1)=-AMP3*PHI0L*RH / PIOMS2 IF(ND.EQ.1) RETURN PHI1C = PHI0C+DEL1C DFI(2)= AMP3*(PHI1C+2D0*PHI0L) / PIOMS2 IF(ND.EQ.2) RETURN PHI1L = PHI0L+DEL1L DFI(3)=-AMP3*(AMP2*PHI1L*RH+6D0*PHI0T/RH) / PIOMS2 IF(ND.EQ.3) RETURN DFI(4)= AMP3*AMP2*(PHI1C+3D0*(DEL1T-DEL1L))/PIOMS2 - 4D0/RH*DFI(3) IF(ND.EQ.4) RETURN END *********************************************************************** SUBROUTINE DFI4DR(ND,RH,AM,ALAM,D0FI,DFI) ** Basic function I_{4} = [-d/dm^{2}+1/Lambda^{2}] I_{2}(m,r) ** and its ND derivatives. ** Each function is scaled by the charged pion mass squared (PIOMS2) IMPLICIT REAL*8 (A-H,O-Z) DIMENSION DFI(*) COMMON/AMSCAL/PIOMS2 DATA SQPI/1.7724538509055160273D0/ RATM = AM/ALAM RATM2= RATM*RATM EXPM = FDEXP(RATM2) XLAM = 0.5D0*RH*ALAM X = AM*RH ERFCM = FDERFC(RATM-XLAM) ERFCP = FDERFC(RATM+XLAM) EPX = 0D0 IF(ERFCP.NE.0) EPX = FDEXP(X) EMX = FDEXP(-X) D0FI = 0.25D0*EXPM*(EMX*ERFCM+EPX*ERFCP)/AM / PIOMS2 DFI(1)=-0.25D0*EXPM*(EMX*ERFCM-EPX*ERFCP) / PIOMS2 IF(ND.EQ.1) RETURN AM2 = AM*AM XLAM2= XLAM*XLAM ELAM = FDEXP(-XLAM2) FORM= AM*ELAM/(2D0*SQPI*RATM * PIOMS2) DFI(2) = D0FI*AM2 - FORM IF(ND.EQ.2) RETURN DFI(3) = DFI(1)*AM2 + FORM*ALAM*XLAM IF(ND.EQ.3) RETURN END *********************************************************************** SUBROUTINE I2NYM(RH,AM,ALAM,PLAM,FI) ** Basic function I_{2} = m phi^{0}_{C} of the Nijmegen model. ** The function is scaled by the charged pion mass squared (PIOMS2) IMPLICIT REAL*8 (A-H,O-Z) COMMON/AMSCAL/PIOMS2 DATA ZLIM/1D-2/, RATLIM/1D0/, SQPI/1.7724538509055160273D0/ AMP =DSQRT(AM*AM+PLAM*PLAM) RATM=AMP/ALAM EXPM=FDEXP(AM*AM/(ALAM*ALAM)) EXPL=FDEXP(-PLAM*PLAM/(ALAM*ALAM)) Z=RH*ALAM IF(Z.GT.ZLIM .OR. RATM.GT.RATLIM) THEN XLAM = 0.5D0*RH*ALAM XLAM2= XLAM*XLAM ELAM = FDEXP(-XLAM2) X = AMP*RH ERFCM = FDERFC(RATM-XLAM) ERFCP = FDERFC(RATM+XLAM) EPX = 0D0 IF(ERFCP.NE.0) EPX = FDEXP(X) EMX = FDEXP(-X) PHI0C=EXPM*(EMX*ERFCM-EPX*ERFCP)/(2D0*X) ELSE A=RATM B=-FDERFC(A)*A*EXPM*SQPI A2=A*A A4=A2*A2 A6=A4*A2 A8=A6*A2 Z2=Z*Z Z4=Z2*Z2 Z6=Z4*Z2 Z8=Z6*Z2 PHI0C = (16*A8*B+16*A8-8*A6+12*A4-30*A2+105)*Z8/5806080D0+ . (8*A6*B+8*A6-4*A4+6*A2-15)*Z6/40320D0+ . (4*A4*B+4*A4-2*A2+3)*Z4/480D0+ . (2*A2*B+2*A2-1)*Z2/12D0+B+1D0 PHI0C = PHI0C/(SQPI*A) * EXPL ENDIF FI = AMP*PHI0C / PIOMS2 RETURN END *********************************************************************** SUBROUTINE I2PHOT(RH,ALAM,PLAM,FI) ** Basic function I_{2} = lim (m-->0) m phi^{0}_{C} ** (for massless photon), scaled by the charged pion mass squared IMPLICIT REAL*8 (A-H,O-Z) COMMON/AMSCAL/PIOMS2 AMP = PLAM RATM = AMP/ALAM XLAM = 0.5D0*RH*ALAM X = AMP*RH ERFCM = FDERFC(RATM-XLAM) ERFCP = FDERFC(RATM+XLAM) EPX = 0D0 IF(ERFCP.NE.0) EPX = FDEXP(X) EMX = FDEXP(-X) PHI0C = (EMX*ERFCM-EPX*ERFCP)/(2D0*RH) FI = PHI0C / PIOMS2 RETURN END *********************************************************************** SUBROUTINE PAIR(R,TYPE,ISO) C* Two-meson-exchange potential form pair terms (seagull vertices). C* Reference : Rijken & Stoks, Phys.Rev.C54 (1996) 2869. C* Mean pion mass everywhere C* C* INPUT: R in fermi C* ----- TYPE 'PP', 'NN', 'NP', 'PN' C* ISO total isospin C* C* OUTPUT: COMMON/VPAIR/V1C,V1S,V1T,V1L, V2C,V2S,V2T,V2L C* ------ C* Pair interactions considered: C* J(PC)=0(++): (\pi\pi)_{0} (\sigma\sigma) (\pi\eta) C* J(PC)=1(--): (\pi\pi)_{1} [g and f couplings] C* J(PC)=1(++): (\pi\rho)_{1} (\pi\sigma) (\pi Pomeron) C* C* IPV = 1: 1/M corrections from vertex functions C* INA = 1: k1.k2/M corrections from nonadiabatic expansion C* = 2: q.(k1-k2)/M corrections from nonadiabatic expansion C* C* Note that the total r-dependent function is divided by PIOMS4 !! C----------------------------------------------------------------------- IMPLICIT REAL*8 (A-H,O-Z) CHARACTER TYPE*2 REAL*8 FINT(8),F2PI(3),I2POM,F2ET(2),F2RO(1),F2ROT(1), . F2SI1(2),F2SI2(2),F2SI(2),F2RO1(1),F2RO2(1) COMMON/VPAIR/ V1C,V1S,V1T,V1L, V2C,V2S,V2T,V2L COMMON/GPAIR/ GPIETA,GPIETP,GSISI, GPISIG,GPIPOM, . GPIPI0,GPIPI1,FPIPI1,GPIRO0,GPIRO1,GPIOM COMMON/MESONM/AMPI,AMETA,AMETP, AMRO,AMOM,AMFI, AMA0,AMEP,AMF0, . AMPIC,AMROC,AMA0C,AWPIC,AWROC,AWA0C, AVSC COMMON/MEANM/ PIM,ROM,A0M COMMON/CUTOFF/ALAM1,ALAMP1,ALAMV1,ALAMS1, . ALAM0,ALAMP0,ALAMV0,ALAMS0 COMMON/CUTTEN/ALAMT1,ALAVT1,ALAMT0,ALAVT0 COMMON/COPLNG/ALPV,THPV,PV1, FPI,FETA,FETP,FPI2,FETA2,FETP2, . ALVE,THV ,GV1, GRO,GOM,GFI, GRO2,GOM2,GFI2, . ALVM, FV1, FRO,FOM,FFI, FRO2,FOM2,FFI2, . ALGS,THGS,GS1, GA0,GEP,GF0, GA02,GEP2,GF02, . GFPRO,GFPOM,GFPFI, GFMRO,GFMOM,GFMFI, . FPIC2,GROC2,FROC2,GFPROC,GFMROC,GA0C2 COMMON/KAPPAV/AKRO,AKOM,AKFI COMMON/YUKEFF/ARO,AMR1,AROC,AMRC1,AWRC1,BRO,AMR2,BROC,AMRC2,AWRC2, . AVSC1,AVSC2, AEPS,AME1,BEPS,AME2 COMMON/POMRON/GPOM,FPOM2,AMPOM,AMPOM2,AMPOM4, . GA2D,FA2D2,AMA2D,AMA2D2,AMA2D4 COMMON/POMFUN/I2POM,DI2POM,DDI2POM COMMON/I2NADI/CSSNA,CPPNA,SOPPNA,CPENA,SOPENA, . SPSNA,TPSNA,SPONA,TPONA COMMON/SCALIN/AMT,AMPV,AMT24 COMMON/AMSCAL/PIOMS2,PIOMS COMMON/HBARC/ HBC DATA SQPI/1.7724538509055160273D0/, TWONPI/0.636619772367581343D0/ DATA NP/12/, IPV/1/, INA/2/ FISO=4D0*ISO-3D0 RH = R/HBC RH2= RH*RH C* General functions: CALL DFI2DR(3,RH,PIM,ALAMP1,0D0,F0PI,F2PI) CALL DFI2DR(2,RH,AME1,ALAMS0,0D0,F0SI1,F2SI1) CALL DFI2DR(2,RH,AME2,ALAMS0,0D0,F0SI2,F2SI2) F0SI =AEPS*F0SI1 +BEPS*F0SI2 F2SI(1)=AEPS*F2SI1(1)+BEPS*F2SI2(1) F2SI(2)=AEPS*F2SI1(2)+BEPS*F2SI2(2) ce CALL DFI2DR(2,RH,AMETA,ALAMP0,0D0,F0ETA,F2ET) CALL DFI2DR(1,RH,AMR1,ALAMV1,0D0,F0RO1,F2RO1) CALL DFI2DR(1,RH,AMR2,ALAMV1,0D0,F0RO2,F2RO2) F0RO =ARO*F0RO1 +BRO*F0RO2 F2RO(1)=ARO*F2RO1(1)+BRO*F2RO2(1) CALL DFI2DR(1,RH,AMR1,ALAVT1,0D0,F0RO1,F2RO1) CALL DFI2DR(1,RH,AMR2,ALAVT1,0D0,F0RO2,F2RO2) F2ROT(1)=ARO*F2RO1(1)+BRO*F2RO2(1) cp XPOM=AMPOM*RH cp XP2=XPOM*XPOM cp I2POM=0D0 cp IF(XP2.LT.170.) I2POM = 4D0*AMPOM/SQPI * (AMPOM/AMT)**2 cp . * FDEXP(-XP2)/PIOMS2 cp DI2POM =-2D0*AMPOM2*RH*I2POM cp DDI2POM=-2D0*AMPOM2*(1D0-2D0*XP2)*I2POM ALAM2= ALAMP1*ALAMP1 PPI0 = (ALAM2/PIOMS2)*ALAMP1/(2D0*SQPI)*FDEXP(-ALAM2*RH2/4D0) CSSNA=F2SI(1)**2 CPPNA=F2PI(2)**2+2D0/RH2*F2PI(1)**2 SOPPNA=2D0/RH2*F2PI(1)**2 ce CPENA=F2PI(2)*F2ET(2)+2D0/RH2*F2PI(1)*F2ET(1) ce SOPENA=2D0/RH2*F2PI(1)*F2ET(1) SPSNA=(F2PI(3)+2D0*(F2PI(2)-F2PI(1)/RH)/RH)*F2SI(1) . -(F2PI(2)*F2SI(2)+2D0*F2PI(1)*F2SI(1)/RH2) TPSNA=(F2PI(3)-(F2PI(2)-F2PI(1)/RH)/RH)*F2SI(1) . -(F2PI(2)*F2SI(2)-F2PI(1)*F2SI(1)/RH2) cp SPONA=(F2PI(3)+2D0*(F2PI(2)-2D0*F2PI(1)/RH)/RH)*DI2POM cp . -F2PI(2)*DDI2POM cp TPONA=(F2PI(3)-(F2PI(2)-2D0*F2PI(1)/RH)/RH)*DI2POM cp . -F2PI(2)*DDI2POM C* Lambda integrals: CALL PAIRINT(RH,1,NP,FINT) B11PI1=PIM*FINT(1)*TWONPI B11PI0=PIM*FINT(2)*TWONPI FINTPI=PIM*FINT(3)*TWONPI B00PI1=PIM*FINT(4)*TWONPI SISINA=PIM*FINT(5)*TWONPI B11SI =PIM*FINT(6)*TWONPI B11RO =PIM*FINT(7)*TWONPI ce B11ETA=PIM*FINT(8)*TWONPI IF(IPV.NE.0 .OR. INA.NE.0) THEN CALL PAIRINT(RH,2,NP,FINT) PI0PVC =PIM*FINT(1)*TWONPI PI0PVL =PIM*FINT(2)*TWONPI PI0NAC =PIM*FINT(3)*TWONPI PI0NAL =PIM*FINT(4)*TWONPI ce PIETPVC=PIM*FINT(5)*TWONPI ce PIETPVL=PIM*FINT(6)*TWONPI ce PIETNAC=PIM*FINT(7)*TWONPI ce PIETNAL=PIM*FINT(8)*TWONPI CALL PAIRINT(RH,3,NP,FINT) PISIPVS=PIM*FINT(1)*TWONPI PISIPVT=PIM*FINT(2)*TWONPI PISINAS=PIM*FINT(3)*TWONPI PISINAT=PIM*FINT(4)*TWONPI cp PIPOPVS=PIM*FINT(5)*TWONPI cp PIPOPVT=PIM*FINT(6)*TWONPI cp PIPONAS=PIM*FINT(7)*TWONPI cp PIPONAT=PIM*FINT(8)*TWONPI ENDIF CALL PAIRINT(RH,4,NP,FINT) B11SIS=PIM*FINT(1)*TWONPI B11SIT=PIM*FINT(2)*TWONPI cp B11POS=PIM*FINT(3)*TWONPI cp B11POT=PIM*FINT(4)*TWONPI C*** One-pair potentials: C*** -------------------- C* Scalar 0++ FF1=GPIPI0*PIOMS*FPI*FPI FF2=GSISI*PIOMS*GEP*GEP*PIOMS2 FF3=FISO*GPIETA*PIOMS*FPI*FETA VC0PP = 6D0*FF1*F2PI(1)*F2PI(1) + 2D0*FF2*F0SI*F0SI ce . + 2D0*FF3*F2PI(1)*F2ET(1) VL0PP = 0D0 IF(IPV.EQ.1) THEN ce VC0PP = VC0PP + 6D0*FF1/AMT*PI0PVC + FF3/AMT*PIETPVC ce VL0PP = 6D0*FF1/AMT*PI0PVL + FF3/AMT*PIETPVL VC0PP = VC0PP + 6D0*FF1/AMT*PI0PVC VL0PP = 6D0*FF1/AMT*PI0PVL ENDIF IF(INA.GE.1) THEN ce VC0PP = VC0PP - (3D0*FF1*PI0NAC+FF3*PIETNAC-FF2*SISINA)/AMT ce IF(INA.EQ.2) VL0PP = VL0PP - (3D0*FF1*PI0NAL+FF3*PIETNAL)/AMT VC0PP = VC0PP - (3D0*FF1*PI0NAC-FF2*SISINA)/AMT IF(INA.EQ.2) VL0PP = VL0PP - 3D0*FF1*PI0NAL/AMT ENDIF C* Vector 1-- FF1=FISO*GPIPI1*FPI*FPI FF2=FISO*4D0/3D0*(GPIPI1+FPIPI1)/AMT*FPI*FPI VC1MM = FF1*4D0*B11PI1 VS1MM =-FF2*F2PI(1)/RH*(F2PI(1)/RH+2D0*F2PI(2)) VT1MM =-FF2*F2PI(1)/RH*(F2PI(1)/RH-F2PI(2)) VL1MM =-FF1*4D0/AMT*F2PI(1)**2/RH2 IF(IPV.EQ.1) THEN VC1MM = VC1MM - FF1/AMT*2D0*(F2PI(2)+2D0*F2PI(1)/RH)*PPI0 VL1MM = VL1MM - FF1/AMT*4D0/RH*F2PI(1)*PPI0 ENDIF IF(INA.GE.1) THEN VC1MM = VC1MM - FF1*2D0/AMT*(F2PI(2)**2+2D0/RH2*F2PI(1)**2) IF(INA.EQ.2) VL1MM = VL1MM - FF1*4D0/AMT*F2PI(1)**2/RH2 ENDIF C* Axial 1++ FF1=FISO/3D0*GPISIG*GEP*PIOMS*FPI FF2=FISO/3D0*GPIPOM*GPOM*PIOMS*FPI VS1PP = FF1*2D0*((F2PI(2)+2D0*F2PI(1)/RH)*F0SI . -F2PI(1)*F2SI(1)) cp . - FF2*2D0*((F2PI(2)+2D0*F2PI(1)/RH)*I2POM-F2PI(1)*DI2POM) VT1PP = FF1*2D0*((F2PI(2)-F2PI(1)/RH)*F0SI . -F2PI(1)*F2SI(1)) cp . - FF2*2D0*((F2PI(2)-F2PI(1)/RH)*I2POM-F2PI(1)*DI2POM) IF(IPV.EQ.1) THEN cp VS1PP = VS1PP + (FF1*PISIPVS-FF2*PIPOPVS)/AMT cp VT1PP = VT1PP + (FF1*PISIPVT-FF2*PIPOPVT)/AMT VS1PP = VS1PP + FF1*PISIPVS/AMT VT1PP = VT1PP + FF1*PISIPVT/AMT ENDIF IF(INA.GE.1) THEN cp VS1PP = VS1PP - (FF1*PISINAS-FF2*PIPONAS)/AMT cp VT1PP = VT1PP - (FF1*PISINAT-FF2*PIPONAT)/AMT VS1PP = VS1PP - FF1*PISINAS/AMT VT1PP = VT1PP - FF1*PISINAT/AMT ENDIF FF3=FISO*2D0/3D0*GPIRO1*PIOMS2*FPI*GRO/AMT VS1PP = VS1PP - FF3*( (F2PI(2)+2D0*F2PI(1)/RH)*F0RO . -2D0*F2PI(1)*(F2RO(1)+AKRO*F2ROT(1)) ) VT1PP = VT1PP - FF3*( (F2PI(2)-F2PI(1)/RH)*F0RO . +F2PI(1)*(F2RO(1)+AKRO*F2ROT(1)) ) V1C = VC0PP+VC1MM V1S = VS1MM+VS1PP V1T = VT1MM+VT1PP V1L = VL0PP+VL1MM C*** Two-pair potentials: C*** -------------------- C* Scalar 0++ WC0PP = -3D0*GPIPI0*GPIPI0*PIOMS2*B11PI0 . -GSISI*GSISI*PIOMS2*B11SI ce . -0.5D0*FISO*GPIETA*GPIETA*PIOMS2*B11ETA C* Vector 1-- WC1MM =-FISO*GPIPI1*GPIPI1*(PPI0*FINTPI-2D0*B00PI1) C* Axial 1++ cp WS1PP = FISO/6D0*(GPISIG*GPISIG*B11SIS - GPIPOM*GPIPOM*B11POS) cp WT1PP = FISO/6D0*(GPISIG*GPISIG*B11SIT - GPIPOM*GPIPOM*B11POT) WS1PP = FISO/6D0*GPISIG*GPISIG*B11SIS WT1PP = FISO/6D0*GPISIG*GPISIG*B11SIT WS1PP = WS1PP - FISO*GPIRO1*GPIRO1*PIOMS2*B11RO V2C = WC0PP+WC1MM V2S = WS1PP V2T = WT1PP V2L = 0D0 RETURN END ************************************************************************ SUBROUTINE PAIRINT(RH,NC,NP,SUM) ** Calculates integrals for pair terms using Gauss integration ** on (0,infinity) to (-1,1) for hyperbolic mapping with NP points IMPLICIT REAL*8 (A-H,O-Z) DIMENSION SUM(8),FUNCW(8),X(20),W(20),NDIM(4) cpce DATA ICAL/0/, NDIM/8,8,8,4/ DATA ICAL/0/, NDIM/7,4,4,2/ SAVE X,W IF(ICAL.EQ.0) THEN CALL GAUSS(NP,X,W,1) DO 1 I=1,NP W(I)=2D0*W(I)/(1D0-X(I))**2 X(I)=(1D0+X(I))/(1D0-X(I)) 1 CONTINUE ICAL=1 ENDIF KDIM=NDIM(NC) DO 10 K=1,KDIM 10 SUM(K)=0D0 DO 11 IY=1,NP CALL FUNCP(FUNCW,X(IY),RH,NC) DO 12 K=1,KDIM 12 SUM(K)=SUM(K)+FUNCW(K)*W(IY) 11 CONTINUE RETURN END ************************************************************************ SUBROUTINE FUNCP(FINT,Y,RH,NC) IMPLICIT REAL*8 (A-H,O-Z) REAL*8 FINT(8),F2PI(2),F3PI(3),I2POM,F2SI(2),F2ET(2), . F1SI1(1),F1SI2(1),F2SI1(2),F2SI2(2) COMMON/CUTOFF/ALAM1,ALAMP1,ALAMV1,ALAMS1, . ALAM0,ALAMP0,ALAMV0,ALAMS0 COMMON/MESONM/AMPI,AMETA,AMETP, AMRO,AMOM,AMFI, AMA0,AMEP,AMF0, . AMPIC,AMROC,AMA0C,AWPIC,AWROC,AWA0C, AVSC COMMON/YUKEFF/ARO,AMR1,AROC,AMRC1,AWRC1,BRO,AMR2,BROC,AMRC2,AWRC2, . AVSC1,AVSC2, AEPS,AME1,BEPS,AME2 COMMON/POMRON/GPOM,FPOM2,AMPOM,AMPOM2,AMPOM4, . GA2D,FA2D2,AMA2D,AMA2D2,AMA2D4 COMMON/POMFUN/I2POM,DI2POM,DDI2POM COMMON/MEANM/ PIM,ROM,A0M COMMON/AMSCAL/PIOMS2,PIOMS COMMON/I2NADI/CSSNA,CPPNA,SOPPNA,CPENA,SOPENA, . SPSNA,TPSNA,SPONA,TPONA RH2=RH*RH PLAM=Y*PIM PLAM2=PLAM*PLAM IF(NC.EQ.1) THEN CALL DFI2DR(2,RH,PIM,ALAMP1,PLAM,F0PI,F2PI) CALL DFI2DR(1,RH,AME1,ALAMS0,PLAM,F0SI1,F1SI1) CALL DFI2DR(1,RH,AME2,ALAMS0,PLAM,F0SI2,F1SI2) F0SI =AEPS*F0SI1 +BEPS*F0SI2 F2SI(1)=AEPS*F1SI1(1)+BEPS*F1SI2(1) CALL I2NYM(RH,AMR1,ALAMV1,PLAM,F0RO1) CALL I2NYM(RH,AMR2,ALAMV1,PLAM,F0RO2) F0RO=ARO*F0RO1+BRO*F0RO2 ce CALL I2NYM(RH,AMETA,ALAMP0,PLAM,F0ETA) ** -(k1.k2)B_{1,1}(pi,pi)_{1} [one pair] FINT(1)=F2PI(1)*F2PI(1) ** B_{1,1}(pi,pi)_{0} [two pair] FINT(2)=F0PI*F0PI ** F(pi)_{1}, B_{0,0}(pi,pi)_{1} [two pair] FINT(3)=F0PI FINT(4)=PLAM2*F0PI*F0PI ** (k_{1}.k_{2}) [ I I - F F ](sigma,sigma) [one pair] FINT(5)=-(CSSNA-F2SI(1)**2)/PLAM2 ** B_{1,1}(sigma,sigma) [two pair] FINT(6)=F0SI*F0SI ** B_{1,1}(pi,rho) [two pair] FINT(7)=F0PI*F0RO ** B_{1,1}(pi,eta) [two pair] ce FINT(8)=F0PI*F0ETA ELSEIF(NC.EQ.2) THEN ** Pseudovector and nonadiabatic for pi-pi_{0}, pi-eta: CALL DFI2DR(2,RH,PIM,ALAMP1,PLAM,F0PI,F2PI) ce CALL DFI2DR(2,RH,AMETA,ALAMP0,PLAM,F0ET,F2ET) ** 1/2(k_{1}^{2}+k_{2}^{2}) B_{1,1}(pi,pi) FINT(1)=-(F2PI(2)+2D0*F2PI(1)/RH)*F0PI ** -i/2(sigma1+sigma2).(k_{1}+k_{2})xq B_{1,1}(pi,pi) FINT(2)=-2D0/RH*F2PI(1)*F0PI ** (k_{1}.k_{2})^{2} [ I I - F F ](pi,pi) FINT(3)=(CPPNA-F2PI(2)**2-2D0*(F2PI(1)/RH)**2)/PLAM2 ** i/2 q.(k1-k2)(sigma1+sigma2).(k1xk2) [ I I - F F ](pi,pi) FINT(4)=(SOPPNA-2D0*(F2PI(1)/RH)**2)/PLAM2 ** (k_{1}^{2}+k_{2}^{2}) B_{1,1}(pi,eta) ce FINT(5)=-(F2PI(2)+2D0*F2PI(1)/RH)*F0ET ce . -F0PI*(F2ET(2)+2D0*F2ET(1)/RH) ** -i(sigma1+sigma2).(k_{1}+k_{2})xq B_{1,1}(pi,eta) ce FINT(6)=-2D0/RH*(F2PI(1)*F0ET+F0PI*F2ET(1)) ** (k_{1}.k_{2})^{2} [ I I - F F ](pi,eta) ce FINT(7)=(CPENA-F2PI(2)*F2ET(2)-2D0/RH2*F2PI(1)*F2ET(1))/PLAM2 ** i/2 q.(k1-k2)(sigma1+sigma2).(k1xk2) [ I I - F F ](pi,eta) ce FINT(8)=(SOPENA-2D0/RH2*F2PI(1)*F2ET(1))/PLAM2 ELSEIF(NC.EQ.3) THEN ** Pseudovector and nonadiabatic for pi-sigma, pi-pomeron: CALL DFI2DR(3,RH,PIM,ALAMP1,PLAM,F0PI,F3PI) CALL DFI2DR(2,RH,AME1,ALAMS0,PLAM,F0SI1,F2SI1) CALL DFI2DR(2,RH,AME2,ALAMS0,PLAM,F0SI2,F2SI2) F2SI(1)=AEPS*F2SI1(1)+BEPS*F2SI2(1) F2SI(2)=AEPS*F2SI1(2)+BEPS*F2SI2(2) ** 3/2[si1.k1si2.k2+si1.k2si2.k1-2si1.k2si2.k2]B_{1,1}(pi,sig) FINT(1)=F0PI*(F2SI(2)+2D0*F2SI(1)/RH)-F3PI(1)*F2SI(1) FINT(2)=F0PI*(F2SI(2)-F2SI(1)/RH)-F3PI(1)*F2SI(1) ** 3/2k1.k2[2si1.k1si2.k1-si1.k1si2.k2-si1.k2si2.k1](II-FF)(pi,sig) SPSF=(F3PI(3)+2D0*(F3PI(2)-F3PI(1)/RH)/RH)*F2SI(1) . -(F3PI(2)*F2SI(2)+2D0*F3PI(1)*F2SI(1)/RH2) TPSF=(F3PI(3)-(F3PI(2)-F3PI(1)/RH)/RH)*F2SI(1) . -(F3PI(2)*F2SI(2)-F3PI(1)*F2SI(1)/RH2) FINT(3)=(SPSNA-SPSF) / PLAM2 FINT(4)=(TPSNA-TPSF) / PLAM2 cp FACPOM=FDEXP(-PLAM2/(4D0*AMPOM2)) ** 3/2[si1.k1si2.k2+si1.k2si2.k1-2si1.k2si2.k2]B_{1,1}(pi,Pom) cp FINT(5)=FACPOM * (F0PI*(DDI2POM+2D0*DI2POM/RH)-F3PI(1)*DI2POM) cp FINT(6)=FACPOM * (F0PI*(DDI2POM-DI2POM/RH)-F3PI(1)*DI2POM) ** 3/2k1.k2[2si1.k1si2.k1-si1.k1si2.k2-si1.k2si2.k1](II-FF)(pi,Pom) cp FINT(7)=(SPONA-((F3PI(3)+2D0*(F3PI(2)-2D0*F3PI(1)/RH)/RH) cp . *DI2POM-F3PI(2)*DDI2POM)*FACPOM) / PLAM2 cp FINT(8)=(TPONA-((F3PI(3)-(F3PI(2)-2D0*F3PI(1)/RH)/RH) cp . *DI2POM-F3PI(2)*DDI2POM)*FACPOM) / PLAM2 ELSEIF(NC.EQ.4) THEN CALL DFI2DR(2,RH,PIM,ALAMP1,PLAM,F0PI,F2PI) CALL DFI2DR(2,RH,AME1,ALAMS0,PLAM,F0SI1,F2SI1) CALL DFI2DR(2,RH,AME2,ALAMS0,PLAM,F0SI2,F2SI2) F0SI =AEPS*F0SI1 +BEPS*F0SI2 F2SI(1)=AEPS*F2SI1(1)+BEPS*F2SI2(1) F2SI(2)=AEPS*F2SI1(2)+BEPS*F2SI2(2) ** 3 sigma1.(k1-k2)sigma2.(k1-k2)B_{1,1}(pi,sigma) FINT(1)=(F2PI(2)+2D0*F2PI(1)/RH)*F0SI . - 2D0*F2PI(1)*F2SI(1) . + F0PI*(F2SI(2)+2D0*F2SI(1)/RH) FINT(2)=(F2PI(2)-F2PI(1)/RH)*F0SI . - 2D0*F2PI(1)*F2SI(1) . + F0PI*(F2SI(2)-F2SI(1)/RH) ** 3 sigma1.(k1-k2)sigma2.(k1-k2)B_{1,1}(pi,Pomeron) cp FACPOM=FDEXP(-PLAM2/(4D0*AMPOM2)) cp FINT(3)=( (F2PI(2)+2D0*F2PI(1)/RH)*I2POM - 2D0*F2PI(1)*DI2POM cp . + F0PI*(DDI2POM+2D0*DI2POM/RH) ) * FACPOM cp FINT(4)=( (F2PI(2)-F2PI(1)/RH)*I2POM - 2D0*F2PI(1)*DI2POM cp . + F0PI*(DDI2POM-DI2POM/RH) ) * FACPOM ENDIF RETURN END ************************************************************************ SUBROUTINE PIMES(R,TYPE,ISO) C* Pion-meson-exchange potentials. C* Reference : Rijken & Stoks, Phys.Rev.C54 (1996) 2851. C* Mean pion mass everywhere C* C* INPUT: R in fermi C* ----- TYPE 'PP', 'NN', 'NP', 'PN' C* ISO total isospin C* C* OUTPUT: COMMON/VPIMES/VCPI,VSPI,VTPI,VLPI, VCPM,VSPM,VTPM,VLPM C* ------ C* Note that the total r-dependent function is divided by PIOMS4 !! C* Standard: leading order Brueckner-Watson (parallel and crossed) C* and Taketani-Machida-Ohnuma (parallel) C* IPV = 1 : 1/M contributions from pseudovector correction C* =-1 : pseudoscalar coupling ==> sign change for INA=3 C* this has to be set also in OM2VK and OM2SCA C* INA = 1 : k1.k2/M contributions from nonadiabatic expansion C* = 2 : q.(k1-k2)/M contributions from nonadiabatic expansion C* = 3 : 1/M off-shell contribution in TMO diagrams C* this has to be set also in OM2VK and OM2SCA C* IM2 = 1 : 1/M^2 contributions from (1+f/g) in vector mesons C* = 2 : 1/M^2 contributions remaining in vector mesons C* = 3 : 1/M^2 contributions in scalar and diffractive mesons C----------------------------------------------------------------------- IMPLICIT REAL*8 (A-H,O-Z) CHARACTER TYPE*2 REAL*8 FINT(11) REAL*8 I2PI,I2ET,I2EP, I2RO,I2OM,I2FI, I2A0,I2SI,I2F0, . I4PI,I4ET,I4EP, I4RO,I4OM,I4FI, I4A0,I4SI,I4F0, . I2POM,I2A2D,I4POM,I4A2D REAL*8 D2PI(4),D2ET(3),D2EP(3), D4PI(3),D4ET(3),D4EP(3), . D2RO(2),D2OM(2),D2FI(2), D4RO(1),D4OM(1),D4FI(1), . D2A0(2),D2SI(1),D2F0(1), D4A0(1),D4SI(1),D4F0(1), . F2RO1(2),F2RO2(2),F2SI1(1),F2SI2(1), . F4RO1(1),F4RO2(1),F4SI1(1),F4SI2(1) REAL*8 DROK(1),DRO2K(2),DROK2(2),DOMK(1),DFIK(1),DOMK2(2),DFIK2(2) common/hoiho/scal COMMON/VPIMES/VCPI,VSPI,VTPI,VLPI, VCPM,VSPM,VTPM,VLPM COMMON/MESONM/AMPI,AMETA,AMETP, AMRO,AMOM,AMFI, AMA0,AMEP,AMF0, . AMPIC,AMROC,AMA0C,AWPIC,AWROC,AWA0C, AVSC COMMON/MEANM/ PIM,ROM,A0M COMMON/COPLNG/ALPV,THPV,PV1, FPI,FETA,FETP,FPI2,FETA2,FETP2, . ALVE,THV ,GV1, GRO,GOM,GFI, GRO2,GOM2,GFI2, . ALVM, FV1, FRO,FOM,FFI, FRO2,FOM2,FFI2, . ALGS,THGS,GS1, GA0,GEP,GF0, GA02,GEP2,GF02, . GFPRO,GFPOM,GFPFI, GFMRO,GFMOM,GFMFI, . FPIC2,GROC2,FROC2,GFPROC,GFMROC,GA0C2 COMMON/KAPPAV/AKRO,AKOM,AKFI COMMON/YUKEFF/ARO,AMR1,AROC,AMRC1,AWRC1,BRO,AMR2,BROC,AMRC2,AWRC2, . AVSC1,AVSC2, AEPS,AME1,BEPS,AME2 COMMON/POMRON/GPOM,FPOM2,AMPOM,AMPOM2,AMPOM4, . GA2D,FA2D2,AMA2D,AMA2D2,AMA2D4 COMMON/SCALIN/AMT,AMPV,AMT24 COMMON/AMSCAL/PIOMS2,PIOMS COMMON/CUTOFF/ALAM1,ALAMP1,ALAMV1,ALAMS1, . ALAM0,ALAMP0,ALAMV0,ALAMS0 COMMON/CUTTEN/ALAMT1,ALAVT1,ALAMT0,ALAVT0 COMMON/IITRMS/PPIIC,PPIIS,PPIIT, . ETIIS,ETIIT, EPIIS,EPIIT, . A0IIS,A0IIT, A2IIS,A2IIT,I2A2D,DI2A2D,DDI2A2D, . RIIEES,RIIEET,RIIEMT,RIIEML,RIIMMC,RIIMMS,RIIMMT, . OIIEMT,OIIMMS,OIIMMT, FIIEMT,FIIMMS,FIIMMT, . ROIIM2S,ROIIM2T,ROIIM2L, ROIIS,ROIIT,ROIIS1,ROIIT1, . OMIIS1,OMIIT1, FIIIS1,FIIIT1, . A0IIM2S,A0IIM2T,A0IIM2L, A2IIM2S,A2IIM2T,A2IIM2L COMMON/HBARC/ HBC DATA SQPI/1.7724538509055160273D0/, TWONPI/0.636619772367581343D0/ DATA SQ2/1.414213562373095D0/ DATA NP/12/, IPV/1/, INA/3/, IM2/3/ ** isospin factors FISO=4D0*ISO-3D0 FPA=3D0-2D0*FISO FCR=3D0+2D0*FISO RH = R/HBC ** I2(m,r) functions CALL DFI2DR(4,RH,PIM, ALAMP1,0D0,I2PI, D2PI) CALL DFI2DR(3,RH,AMETA,ALAMP0,0D0,I2ET, D2ET) CALL DFI2DR(3,RH,AMETP,ALAMP0,0D0,I2EP, D2EP) CALL DFI2DR(2,RH,A0M, ALAMS1,0D0,I2A0, D2A0) CALL DFI2DR(1,RH,AME1, ALAMS0,0D0,F0SI1,F2SI1) CALL DFI2DR(1,RH,AME2, ALAMS0,0D0,F0SI2,F2SI2) CALL DFI2DR(1,RH,AMF0, ALAMS0,0D0,I2F0, D2F0) CALL DFI2DR(2,RH,AMR1, ALAMV1,0D0,F0RO1,F2RO1) CALL DFI2DR(2,RH,AMR2, ALAMV1,0D0,F0RO2,F2RO2) CALL DFI2DR(2,RH,AMOM, ALAMV0,0D0,I2OM, D2OM) CALL DFI2DR(2,RH,AMFI, ALAMV0,0D0,I2FI, D2FI) I2SI =AEPS*F0SI1 +BEPS*F0SI2 D2SI(1)=AEPS*F2SI1(1)+BEPS*F2SI2(1) I2RO =ARO*F0RO1 +BRO*F0RO2 D2RO(1)=ARO*F2RO1(1)+BRO*F2RO2(1) D2RO(2)=ARO*F2RO1(2)+BRO*F2RO2(2) XPOM=AMPOM*RH XP2=XPOM*XPOM I2POM=0D0 IF(XP2.LT.170.) I2POM=FDEXP(-XP2) DI2POM =-2D0*AMPOM2*RH*I2POM DDI2POM=-2D0*AMPOM2*(1D0-2D0*XP2)*I2POM cp XA2=AMA2D*RH cp XP2=XA2*XA2 cp I2A2D=0D0 cp IF(XP2.LT.170.) I2A2D=FDEXP(-XP2) cp DI2A2D =-2D0*AMA2D2*RH*I2A2D cp DDI2A2D=-2D0*AMA2D2*(1D0-2D0*XP2)*I2A2D IF(IPV.NE.0 .OR. INA.NE.0) THEN PIM2 =PIM*PIM AMET2=AMETA*AMETA AMEP2=AMETP*AMETP ** I0'(m,r) functions for 1/M corrections in pseudoscalars XLAM =0.5D0*ALAMP1*RH XLAM2=XLAM*XLAM ALAM2=ALAMP1*ALAMP1 DI0PI=-(ALAM2/PIOMS2)*ALAM2*XLAM/(2D0*SQPI)*FDEXP(-XLAM2) XLAM =0.5D0*ALAMP0*RH XLAM2=XLAM*XLAM ALAM2=ALAMP0*ALAMP0 DI0PS=-(ALAM2/PIOMS2)*ALAM2*XLAM/(2D0*SQPI)*FDEXP(-XLAM2) ENDIF IF(INA.NE.0) THEN ** I4(m,r) functions for nonadiabatic corrections CALL DFI4DR(3,RH,PIM, ALAMP1,I4PI, D4PI) CALL DFI4DR(3,RH,AMETA,ALAMP0,I4ET, D4ET) CALL DFI4DR(3,RH,AMETP,ALAMP0,I4EP, D4EP) CALL DFI4DR(1,RH,A0M, ALAMS1,I4A0, D4A0) CALL DFI4DR(1,RH,AME1, ALAMS0,F0SI1,F4SI1) CALL DFI4DR(1,RH,AME2, ALAMS0,F0SI2,F4SI2) CALL DFI4DR(1,RH,AMF0, ALAMS0,I4F0, D4F0) CALL DFI4DR(1,RH,AMR1, ALAMV1,F0RO1,F4RO1) CALL DFI4DR(1,RH,AMR2, ALAMV1,F0RO2,F4RO2) CALL DFI4DR(1,RH,AMOM, ALAMV0,I4OM, D4OM) CALL DFI4DR(1,RH,AMFI, ALAMV0,I4FI, D4FI) I4SI =AEPS*F0SI1 +BEPS*F0SI2 D4SI(1)=AEPS*F4SI1(1)+BEPS*F4SI2(1) I4RO =ARO*F0RO1 +BRO*F0RO2 D4RO(1)=ARO*F4RO1(1)+BRO*F4RO2(1) I4POM =I2POM /(4D0*AMPOM2) DI4POM=DI2POM/(4D0*AMPOM2) cp I4A2D =I2A2D /(4D0*AMA2D2) cp DI4A2D=DI2A2D/(4D0*AMA2D2) IF(INA.GE.3) THEN ** Common combination of pion function PIS=D2PI(2)+2D0*D2PI(1)/RH PIT=D2PI(2)-D2PI(1)/RH ENDIF ENDIF CALL OPSEU(RH,D2PI(1),D2PI(2),D2PI(1),D2PI(2),PPIIC,PPIIS,PPIIT) CALL TMEINT(RH,1,NP,FINT) BPIPIC = PIM*FINT(1)*TWONPI BPIPIS = PIM*FINT(2)*TWONPI BPIPIT = PIM*FINT(3)*TWONPI CALL OPSEU(RH,D2PI(1),D2PI(2),D2ET(1),D2ET(2),ETIIC,ETIIS,ETIIT) CALL OPSEU(RH,D2PI(1),D2PI(2),D2EP(1),D2EP(2),EPIIC,EPIIS,EPIIT) CALL OSCL(RH,D2PI(1),D2PI(2),I2A0,A0IIS,A0IIT) CALL OVC1(RH,D2PI(1),D2PI(2),I2RO,RIIEES,RIIEET) cp CALL OSCL(RH,D2PI(1),D2PI(2),I2A2D,A2IIS,A2IIT) CALL TMEINT(RH,2,NP,FINT) BPIETS = PIM*FINT(1) *TWONPI BPIETT = PIM*FINT(2) *TWONPI BPIEPS = PIM*FINT(3) *TWONPI BPIEPT = PIM*FINT(4) *TWONPI BPIA0S = PIM*FINT(5) *TWONPI BPIA0T = PIM*FINT(6) *TWONPI BROEES = PIM*FINT(7) *TWONPI BROEET = PIM*FINT(8) *TWONPI cp BPIA2S = PIM*FINT(9) *TWONPI cp BPIA2T = PIM*FINT(10)*TWONPI IF(IM2.GE.1) THEN ** Vector-meson functions for 1+f/g and (1+f/g)**2 CALL SPECV(1,RH,0D0,D2RO,DROK,DROK2,AKRO) CALL SPECV(2,RH,0D0,D2OM,DOMK,DOMK2,AKOM) CALL SPECV(3,RH,0D0,D2FI,DFIK,DFIK2,AKFI) ** Operators for 1/M^2 corrections for (1+f/g) and (1+f/g)**2 CALL OVC2(RH,D2PI(1),D2PI(2),DROK(1),DROK2(1),DROK2(2), . RIIEMT,RIIMMS,RIIMMT) CALL OVC3(RH,D2PI(1),D2PI(2),DROK(1),DROK2(1),DROK2(2), . RIIEML,RIIMMC) CALL OVC2(RH,D2PI(1),D2PI(2),DOMK(1),DOMK2(1),DOMK2(2), . OIIEMT,OIIMMS,OIIMMT) CALL OVC2(RH,D2PI(1),D2PI(2),DFIK(1),DFIK2(1),DFIK2(2), . FIIEMT,FIIMMS,FIIMMT) CALL TMEINT(RH,3,NP,FINT) BROEMT = PIM*FINT(1) *TWONPI/AMT24 BROEML = PIM*FINT(2) *TWONPI/AMT24 BROMMC = PIM*FINT(3) *TWONPI/AMT24 BROMMS = PIM*FINT(4) *TWONPI/AMT24 BROMMT = PIM*FINT(5) *TWONPI/AMT24 BOMEMT = PIM*FINT(6) *TWONPI/AMT24 BOMMMS = PIM*FINT(7) *TWONPI/AMT24 BOMMMT = PIM*FINT(8) *TWONPI/AMT24 BFIEMT = PIM*FINT(9) *TWONPI/AMT24 BFIMMS = PIM*FINT(10)*TWONPI/AMT24 BFIMMT = PIM*FINT(11)*TWONPI/AMT24 ENDIF IF(IM2.GE.2) THEN ** Vector-meson functions for 1+2f/g CALL SPECVR(RH,0D0,D2RO,DRO2K,AKRO) ** Operators for remaining 1/M^2 corrections for vector mesons CALL OM2VEC(RH,D2PI(1),D2PI(2),DRO2K(1),DRO2K(2), . ROIIM2S,ROIIM2T,ROIIM2L) CALL OM2VK(RH,D2PI(1),D2PI(2),D2PI(3),D2PI(4), . I2RO,D2RO(1),D2RO(2), ROIIS,ROIIT) CALL OM2VK1(RH,D2PI(1),D2PI(2),D2PI(3),D2PI(4),I2RO, . ROIIS1,ROIIT1) CALL OM2VK1(RH,D2PI(1),D2PI(2),D2PI(3),D2PI(4),I2OM, . OMIIS1,OMIIT1) CALL OM2VK1(RH,D2PI(1),D2PI(2),D2PI(3),D2PI(4),I2FI, . FIIIS1,FIIIT1) CALL TMEINT(RH,4,NP,FINT) BROM2S = PIM*FINT(1) *TWONPI/AMT24 BROM2T = PIM*FINT(2) *TWONPI/AMT24 BROM2L = PIM*FINT(3) *TWONPI/AMT24 BROSK = PIM*FINT(4) *TWONPI/AMT24 BROTK = PIM*FINT(5) *TWONPI/AMT24 BROSK1 = PIM*FINT(6) *TWONPI/AMT24 BROTK1 = PIM*FINT(7) *TWONPI/AMT24 BOMSK1 = PIM*FINT(8) *TWONPI/AMT24 BOMTK1 = PIM*FINT(9) *TWONPI/AMT24 BFISK1 = PIM*FINT(10)*TWONPI/AMT24 BFITK1 = PIM*FINT(11)*TWONPI/AMT24 ENDIF IF(IM2.EQ.3) THEN ** Operators for 1/M^2 corrections for scalar mesons CALL OM2SCA(RH,D2PI(1),D2PI(2),D2PI(3),D2A0(1),D2A0(2), . A0IIM2S,A0IIM2T,A0IIM2L) cp CALL OM2SCA(RH,D2PI(1),D2PI(2),D2PI(3),DI2A2D,DDI2A2D, cp . A2IIM2S,A2IIM2T,A2IIM2L) CALL TMEINT(RH,5,NP,FINT) BA0M2S = PIM*FINT(1)*TWONPI/AMT24 BA0M2T = PIM*FINT(2)*TWONPI/AMT24 BA0M2L = PIM*FINT(3)*TWONPI/AMT24 cp BA2M2S = PIM*FINT(4)*TWONPI/AMT24 cp BA2M2T = PIM*FINT(5)*TWONPI/AMT24 cp BA2M2L = PIM*FINT(6)*TWONPI/AMT24 ENDIF C* Two-pion exchange: FPI4=FPI*FPI*FPI*FPI VCPI = FPI4*0.5D0*(FPA-FCR)*BPIPIC VSPI =-FPI4*0.5D0*(FPA+FCR)*BPIPIS VTPI =-FPI4*0.5D0*(FPA+FCR)*BPIPIT VLPI = 0D0 IF(IPV.EQ.1) THEN VCPI = VCPI + FCR*FPI4*(D2PI(1)*(PIM2*D2PI(1)-DI0PI))/AMT VLPI = VLPI + FCR*FPI4*(D2PI(2)+D2PI(1)/RH)*D2PI(1)*2D0/RH/AMT ENDIF IF(INA.NE.0) THEN CALL ONAPS(RH,D2PI(1),D2PI(2),D2PI(3),D4PI(1),D4PI(2),D4PI(3), . CPIPI,SPIPI,TPIPI,OPIPI) VCPI = VCPI + FPI4/AMT * ( 0.5D0*FPA-FCR)*CPIPI VSPI = VSPI + FPI4/AMT * (-0.5D0*FPA-FCR)*SPIPI VTPI = VTPI + FPI4/AMT * (-0.5D0*FPA-FCR)*TPIPI IF(INA.GE.2) VLPI = VLPI + FPI4/AMT * 0.5D0*FPA*OPIPI IF(INA.GE.3) THEN ** extra part form iterated OPE diagrams (i.e., TMO diagrams) VCPI = VCPI + IPV*FPA*FPI4 . *(D2PI(1)*(PIM2*D2PI(1)-DI0PI))/2D0/AMT VLPI = VLPI + IPV*FPA*FPI4*(D2PI(2)-D2PI(1)/RH)*D2PI(1)/RH/AMT ENDIF ENDIF VCPM = 0D0 VSPM = 0D0 VTPM = 0D0 VLPM = 0D0 C* Pseudoscalar (eta): FF=FPI*FPI*FETA*FETA * FISO VSPM = VSPM - FF*2D0*BPIETS VTPM = VTPM - FF*2D0*BPIETT IF(IPV.EQ.1) THEN VCPM = VCPM + FF*( (PIM2+AMET2)*D2PI(1)*D2ET(1) . -DI0PI*D2ET(1)-D2PI(1)*DI0PS )/AMT VLPM = VLPM + FF*( (D2PI(2)+D2PI(1)/RH)*D2ET(1) . +(D2ET(2)+D2ET(1)/RH)*D2PI(1) )*2D0/RH/AMT ENDIF IF(INA.NE.0) THEN CALL ONAPS(RH,D2PI(1),D2PI(2),D2PI(3),D4ET(1),D4ET(2),D4ET(3), . CPIET24,SPIET24,TPIET24,OPIET24) CALL ONAPS(RH,D4PI(1),D4PI(2),D4PI(3),D2ET(1),D2ET(2),D2ET(3), . CPIET42,SPIET42,TPIET42,OPIET42) VCPM = VCPM - FF/AMT * (CPIET24+CPIET42)*0.5D0 VSPM = VSPM - FF/AMT * (SPIET24+SPIET42)*1.5D0 VTPM = VTPM - FF/AMT * (TPIET24+TPIET42)*1.5D0 IF(INA.GE.2) VLPM = VLPM + FF/AMT * (OPIET24+OPIET42)*0.5D0 IF(INA.GE.3) THEN ** extra part form iterated diagrams (i.e., TMO diagrams) VCPM = VCPM + IPV*FF*( (PIM2+AMET2)*D2PI(1)*D2ET(1) . -DI0PI*D2ET(1)-D2PI(1)*DI0PS )/2D0/AMT VLPM = VLPM + IPV*FF*( (D2PI(2)-D2PI(1)/RH)*D2ET(1) . +(D2ET(2)-D2ET(1)/RH)*D2PI(1) )/RH/AMT ENDIF ENDIF C* Pseudoscalar (eta'): FF=FPI*FPI*FETP*FETP * FISO VSPM = VSPM - FF*2D0*BPIEPS VTPM = VTPM - FF*2D0*BPIEPT IF(IPV.EQ.1) THEN VCPM = VCPM + FF*( (PIM2+AMEP2)*D2PI(1)*D2EP(1) . -DI0PI*D2EP(1)-D2PI(1)*DI0PS )/AMT VLPM = VLPM + FF*( (D2PI(2)+D2PI(1)/RH)*D2EP(1) . +(D2EP(2)+D2EP(1)/RH)*D2PI(1) )*2D0/RH/AMT ENDIF IF(INA.NE.0) THEN CALL ONAPS(RH,D2PI(1),D2PI(2),D2PI(3),D4EP(1),D4EP(2),D4EP(3), . CPIEP24,SPIEP24,TPIEP24,OPIEP24) CALL ONAPS(RH,D4PI(1),D4PI(2),D4PI(3),D2EP(1),D2EP(2),D2EP(3), . CPIEP42,SPIEP42,TPIEP42,OPIEP42) VCPM = VCPM - FF/AMT * (CPIEP24+CPIEP42)*0.5D0 VSPM = VSPM - FF/AMT * (SPIEP24+SPIEP42)*1.5D0 VTPM = VTPM - FF/AMT * (TPIEP24+TPIEP42)*1.5D0 IF(INA.GE.2) VLPM = VLPM + FF/AMT * (OPIEP24+OPIEP42)*0.5D0 IF(INA.GE.3) THEN ** extra part form iterated diagrams (i.e., TMO diagrams) VCPM = VCPM + IPV*FF*( (PIM2+AMEP2)*D2PI(1)*D2EP(1) . -DI0PI*D2EP(1)-D2PI(1)*DI0PS )/2D0/AMT VLPM = VLPM + IPV*FF*( (D2PI(2)-D2PI(1)/RH)*D2EP(1) . +(D2EP(2)-D2EP(1)/RH)*D2PI(1) )/RH/AMT ENDIF ENDIF C* Scalar (a0): FF=FPI*FPI*GA0*GA0*PIOMS2 VSPM = VSPM + FF*(FPA-FCR)*BPIA0S VTPM = VTPM + FF*(FPA-FCR)*BPIA0T IF(IPV.EQ.1) THEN VSPM = VSPM - FF*FCR*D2PI(1)*D2A0(1) / (3D0*AMT) VTPM = VTPM - FF*FCR*D2PI(1)*D2A0(1) / (3D0*AMT) ENDIF IF(INA.NE.0) THEN CALL ONAVS(RH,D2PI(1),D2PI(2),D2PI(3),D4A0(1),SPIA024,TPIA024) CALL ONAVS(RH,D4PI(1),D4PI(2),D4PI(3),D2A0(1),SPIA042,TPIA042) VSPM = VSPM + FF/AMT * (0.5D0*FPA-FCR)*(SPIA024+SPIA042) VTPM = VTPM + FF/AMT * (0.5D0*FPA-FCR)*(TPIA024+TPIA042) IF(INA.GE.3) THEN VSPM = VSPM + FF*FPA*(PIS*I2A0-IPV*D2PI(1)*D2A0(1))/(6D0*AMT) VTPM = VTPM + FF*FPA*(PIT*I2A0-IPV*D2PI(1)*D2A0(1))/(6D0*AMT) ENDIF ENDIF IF(IM2.GE.3) THEN VSPM = VSPM + FF*(FPA-FCR)*BA0M2S VTPM = VTPM + FF*(FPA-FCR)*BA0M2T VLPM = VLPM + FF*(FPA-FCR)*BA0M2L ENDIF C* Scalar (epsilon): FF=FPI*FPI*GEP*GEP*PIOMS2 * FISO IF(IPV.EQ.1) THEN VSPM = VSPM - FF*D2PI(1)*D2SI(1) / (3D0*AMT) VTPM = VTPM - FF*D2PI(1)*D2SI(1) / (3D0*AMT) ENDIF IF(INA.NE.0) THEN CALL ONAVS(RH,D2PI(1),D2PI(2),D2PI(3),D4SI(1),SPISI24,TPISI24) CALL ONAVS(RH,D4PI(1),D4PI(2),D4PI(3),D2SI(1),SPISI42,TPISI42) VSPM = VSPM - FF/AMT * (SPISI24+SPISI42)/2D0 VTPM = VTPM - FF/AMT * (TPISI24+TPISI42)/2D0 IF(INA.GE.3) THEN VSPM = VSPM + FF*(PIS*I2SI-IPV*D2PI(1)*D2SI(1))/(6D0*AMT) VTPM = VTPM + FF*(PIT*I2SI-IPV*D2PI(1)*D2SI(1))/(6D0*AMT) ENDIF ENDIF C* Scalar (f0): FF=FPI*FPI*GF0*GF0*PIOMS2 * FISO IF(IPV.EQ.1) THEN VSPM = VSPM - FF*D2PI(1)*D2F0(1) / (3D0*AMT) VTPM = VTPM - FF*D2PI(1)*D2F0(1) / (3D0*AMT) ENDIF IF(INA.NE.0) THEN CALL ONAVS(RH,D2PI(1),D2PI(2),D2PI(3),D4F0(1),SPIF024,TPIF024) CALL ONAVS(RH,D4PI(1),D4PI(2),D4PI(3),D2F0(1),SPIF042,TPIF042) VSPM = VSPM - FF/AMT * (SPIF024+SPIF042)/2D0 VTPM = VTPM - FF/AMT * (TPIF024+TPIF042)/2D0 IF(INA.GE.3) THEN VSPM = VSPM + FF*(PIS*I2F0-IPV*D2PI(1)*D2F0(1))/(6D0*AMT) VTPM = VTPM + FF*(PIT*I2F0-IPV*D2PI(1)*D2F0(1))/(6D0*AMT) ENDIF ENDIF C* Vector (rho): FF=FPI*FPI*GRO*GRO*PIOMS2 VSPM = VSPM + FF * (FPA-FCR)*BROEES VTPM = VTPM + FF * (FPA-FCR)*BROEET IF(IPV.EQ.1) THEN VSPM = VSPM + FF*FCR*D2PI(1)*D2RO(1) / (3D0*AMT) VTPM = VTPM + FF*FCR*D2PI(1)*D2RO(1) / (3D0*AMT) ENDIF IF(INA.NE.0) THEN CALL ONAVS(RH,D2PI(1),D2PI(2),D2PI(3),D4RO(1),SPIRO24,TPIRO24) CALL ONAVS(RH,D4PI(1),D4PI(2),D4PI(3),D2RO(1),SPIRO42,TPIRO42) VSPM = VSPM - FF/AMT * (0.5D0*FPA-FCR)*(SPIRO24+SPIRO42) VTPM = VTPM - FF/AMT * (0.5D0*FPA-FCR)*(TPIRO24+TPIRO42) IF(INA.GE.3) THEN VSPM = VSPM - FF*FPA*(PIS*I2RO-IPV*D2PI(1)*D2RO(1))/(6D0*AMT) VTPM = VTPM - FF*FPA*(PIT*I2RO-IPV*D2PI(1)*D2RO(1))/(6D0*AMT) ENDIF ENDIF IF(IM2.GE.1) THEN VCPM = VCPM + FF * (FPA-FCR)*BROMMC VSPM = VSPM + FF * (FPA+FCR)*BROMMS VTPM = VTPM + FF * (FPA+FCR)*(BROEMT+BROMMT) VLPM = VLPM + FF * (FPA-FCR)*BROEML IF(IM2.GE.2) THEN VSPM = VSPM + FF*((FPA-FCR)*(BROSK+BROM2S)-(FPA+FCR)*BROSK1) VTPM = VTPM + FF*((FPA-FCR)*(BROTK+BROM2T)-(FPA+FCR)*BROTK1) VLPM = VLPM + FF*(FPA-FCR)*BROM2L ENDIF ENDIF C* Vector (omega): FF=FPI*FPI*GOM*GOM*PIOMS2 * FISO IF(IPV.EQ.1) THEN VSPM = VSPM + FF*D2PI(1)*D2OM(1) / (3D0*AMT) VTPM = VTPM + FF*D2PI(1)*D2OM(1) / (3D0*AMT) ENDIF IF(INA.NE.0) THEN CALL ONAVS(RH,D2PI(1),D2PI(2),D2PI(3),D4OM(1),SPIOM24,TPIOM24) CALL ONAVS(RH,D4PI(1),D4PI(2),D4PI(3),D2OM(1),SPIOM42,TPIOM42) VSPM = VSPM + FF/AMT * (SPIOM24+SPIOM42)/2D0 VTPM = VTPM + FF/AMT * (TPIOM24+TPIOM42)/2D0 IF(INA.GE.3) THEN VSPM = VSPM - FF*(PIS*I2OM-IPV*D2PI(1)*D2OM(1))/(6D0*AMT) VTPM = VTPM - FF*(PIT*I2OM-IPV*D2PI(1)*D2OM(1))/(6D0*AMT) ENDIF ENDIF IF(IM2.GE.1) THEN VSPM = VSPM + FF * 2D0*BOMMMS VTPM = VTPM + FF * 2D0*(BOMEMT+BOMMMT) IF(IM2.GE.2) THEN VSPM = VSPM - FF * 2D0*BOMSK1 VTPM = VTPM - FF * 2D0*BOMTK1 ENDIF ENDIF C* Vector (phi): FF=FPI*FPI*GFI*GFI*PIOMS2 * FISO IF(IPV.EQ.1) THEN VSPM = VSPM + FF*D2PI(1)*D2FI(1) / (3D0*AMT) VTPM = VTPM + FF*D2PI(1)*D2FI(1) / (3D0*AMT) ENDIF IF(INA.NE.0) THEN CALL ONAVS(RH,D2PI(1),D2PI(2),D2PI(3),D4FI(1),SPIFI24,TPIFI24) CALL ONAVS(RH,D4PI(1),D4PI(2),D4PI(3),D2FI(1),SPIFI42,TPIFI42) VSPM = VSPM + FF/AMT * (SPIFI24+SPIFI42)/2D0 VTPM = VTPM + FF/AMT * (TPIFI24+TPIFI42)/2D0 IF(INA.GE.3) THEN VSPM = VSPM - FF*(PIS*I2FI-IPV*D2PI(1)*D2FI(1))/(6D0*AMT) VTPM = VTPM - FF*(PIT*I2FI-IPV*D2PI(1)*D2FI(1))/(6D0*AMT) ENDIF ENDIF IF(IM2.GE.1) THEN VSPM = VSPM + FF * 2D0*BFIMMS VTPM = VTPM + FF * 2D0*(BFIEMT+BFIMMT) IF(IM2.GE.2) THEN VSPM = VSPM - FF * 2D0*BFISK1 VTPM = VTPM - FF * 2D0*BFITK1 ENDIF ENDIF C* Pi-Pomeron (isosinglet): FF=FPI*FPI*FPOM2 * FISO IF(IPV.EQ.1) THEN VSPM = VSPM + FF*D2PI(1)*DI2POM / (3D0*AMT) VTPM = VTPM + FF*D2PI(1)*DI2POM / (3D0*AMT) ENDIF IF(INA.NE.0) THEN CALL ONAVS(RH,D2PI(1),D2PI(2),D2PI(3),DI4POM,SPIPO24,TPIPO24) CALL ONAVS(RH,D4PI(1),D4PI(2),D4PI(3),DI2POM,SPIPO42,TPIPO42) VSPM = VSPM + FF/AMT * (SPIPO24+SPIPO42)/2D0 VTPM = VTPM + FF/AMT * (TPIPO24+TPIPO42)/2D0 IF(INA.GE.3) THEN VSPM = VSPM - FF*(PIS*I2POM-IPV*D2PI(1)*DI2POM)/(6D0*AMT) VTPM = VTPM - FF*(PIT*I2POM-IPV*D2PI(1)*DI2POM)/(6D0*AMT) ENDIF ENDIF C* Pi-diffractive (a2): cp FF=FPI*FPI*FA2D2 cp VSPM = VSPM - FF*(FPA-FCR)*BPIA2S cp VTPM = VTPM - FF*(FPA-FCR)*BPIA2T cp IF(IPV.EQ.1) THEN cp VSPM = VSPM + FF*FCR*D2PI(1)*DI2A2D / (3D0*AMT) cp VTPM = VTPM + FF*FCR*D2PI(1)*DI2A2D / (3D0*AMT) cp ENDIF cp IF(INA.NE.0) THEN cp CALL ONAVS(RH,D2PI(1),D2PI(2),D2PI(3),DI4A2D,SPIA224,TPIA224) cp CALL ONAVS(RH,D4PI(1),D4PI(2),D4PI(3),DI2A2D,SPIA242,TPIA242) cp VSPM = VSPM - FF/AMT * (0.5D0*FPA-FCR)*(SPIA224+SPIA242) cp VTPM = VTPM - FF/AMT * (0.5D0*FPA-FCR)*(TPIA224+TPIA242) cp IF(INA.GE.3) THEN cp VSPM = VSPM - FF*FPA*(PIS*I2A2D-IPV*D2PI(1)*DI2A2D)/(6D0*AMT) cp VTPM = VTPM - FF*FPA*(PIT*I2A2D-IPV*D2PI(1)*DI2A2D)/(6D0*AMT) cp ENDIF cp ENDIF cp IF(IM2.GE.3) THEN cp VSPM = VSPM - FF*(FPA-FCR)*BA2M2S cp VTPM = VTPM - FF*(FPA-FCR)*BA2M2T cp VLPM = VLPM - FF*(FPA-FCR)*BA2M2L cp ENDIF C* Scalar-scalar exchange (epsilon, pomeron, and diffractive): *** Use: 2/pi int dl/l**2 (1-exp[-l**2/2ma**2])=sqrt(2/pi)/ma cp FF=0.5D0*FA2D2*FA2D2*(FPA-FCR)*(SQ2/SQPI)/AMA2D cp VCPM = VCPM + FF*I2A2D*I2A2D cp IF(IM2.GE.3) THEN cp VCPM = VCPM - FF*2D0*(DDI2A2D+2D0/RH*DI2A2D)*I2A2D / AMT24 cp VLPM = VLPM - FF*4D0*I2A2D*DI2A2D/RH / AMT24 cp ENDIF cp IF(INA.NE.0) THEN cp FFSS=(GEP*GEP*PIOMS2)**2 cp FFSP=-GEP*GEP*PIOMS2*FPOM2 cp FFSA=-GEP*GEP*PIOMS2*FA2D2*FISO cp FFPP= FPOM2*FPOM2 cp FFPA= FPOM2*FA2D2*FISO cp FFAA= FA2D2*FA2D2*(FCR-0.5D0*FPA) cp VCPM = VCPM + ( FFSS*D2SI(1)*D4SI(1) cp . +FFSP*(D2SI(1)*DI4POM+D4SI(1)*DI2POM) cp . +FFSA*(D2SI(1)*DI4A2D+D4SI(1)*DI2A2D) cp . +FFPP*DI2POM*DI4POM cp . +FFPA*(DI2POM*DI4A2D+DI4POM*DI2A2D) cp . +FFAA*DI2A2D*DI4A2D*2D0 ) / (2D0*AMT) cp ENDIF RETURN END ************************************************************************ SUBROUTINE TMEINT(RH,NC,NP,SUM) ** Calculates integrals for pi-meson terms using Gauss integration ** on (0,infinity) to (-1,1) for hyperbolic mapping with NP points IMPLICIT REAL*8 (A-H,O-Z) DIMENSION SUM(11),FUNCW(11),X(20),W(20),NDIM(5) cp DATA ICAL/0/, NDIM/3,10,11,11,6/ DATA ICAL/0/, NDIM/3,8,11,11,3/ SAVE X,W IF(ICAL.EQ.0) THEN CALL GAUSS(NP,X,W,1) DO 1 I=1,NP W(I)=2D0*W(I)/(1D0-X(I))**2 X(I)=(1D0+X(I))/(1D0-X(I)) 1 CONTINUE ICAL=1 ENDIF KDIM=NDIM(NC) DO 10 K=1,KDIM 10 SUM(K)=0D0 DO 11 IY=1,NP CALL FUNCT(FUNCW,X(IY),RH,NC) DO 12 K=1,KDIM 12 SUM(K)=SUM(K)+FUNCW(K)*W(IY) 11 CONTINUE RETURN END ************************************************************************ SUBROUTINE FUNCT(FINT,Y,RH,NC) IMPLICIT REAL*8 (A-H,O-Z) REAL*8 FINT(11), F2PI(2),F2PI3(3),F2PI4(4), . F2ET(2),F2EP(2),F2A0(2), F2RO(2),F2OM(2),F2FI(2), . F2RO1(2),F2RO2(2), I2A2D REAL*8 FROK(1),FRO2K(2),FROK2(2),FOMK(1),FFIK(1),FOMK2(2),FFIK2(2) COMMON/MESONM/AMPI,AMETA,AMETP, AMRO,AMOM,AMFI, AMA0,AMEP,AMF0, . AMPIC,AMROC,AMA0C,AWPIC,AWROC,AWA0C, AVSC COMMON/KAPPAV/AKRO,AKOM,AKFI COMMON/YUKEFF/ARO,AMR1,AROC,AMRC1,AWRC1,BRO,AMR2,BROC,AMRC2,AWRC2, . AVSC1,AVSC2, AEPS,AME1,BEPS,AME2 COMMON/CUTOFF/ALAM1,ALAMP1,ALAMV1,ALAMS1, . ALAM0,ALAMP0,ALAMV0,ALAMS0 COMMON/MEANM/ PIM,ROM,A0M COMMON/POMRON/GPOM,FPOM2,AMPOM,AMPOM2,AMPOM4, . GA2D,FA2D2,AMA2D,AMA2D2,AMA2D4 COMMON/IITRMS/PPIIC,PPIIS,PPIIT, . ETIIS,ETIIT, EPIIS,EPIIT, . A0IIS,A0IIT, A2IIS,A2IIT,I2A2D,DI2A2D,DDI2A2D, . RIIEES,RIIEET,RIIEMT,RIIEML,RIIMMC,RIIMMS,RIIMMT, . OIIEMT,OIIMMS,OIIMMT, FIIEMT,FIIMMS,FIIMMT, . ROIIM2S,ROIIM2T,ROIIM2L, ROIIS,ROIIT,ROIIS1,ROIIT1, . OMIIS1,OMIIT1, FIIIS1,FIIIT1, . A0IIM2S,A0IIM2T,A0IIM2L, A2IIM2S,A2IIM2T,A2IIM2L PLAM=Y*PIM PLAM2=PLAM*PLAM IF(NC.EQ.1) THEN *** FINT = dlam/lam**2 [(I2 operator I2) - (F operator F)](pi,pi) CALL DFI2DR(2,RH,PIM, ALAMP1,PLAM,F0PI, F2PI) CALL OPSEU(RH,F2PI(1),F2PI(2),F2PI(1),F2PI(2),PIFFC,PIFFS,PIFFT) FINT(1) = (PPIIC-PIFFC) / PLAM2 FINT(2) = (PPIIS-PIFFS) / PLAM2 FINT(3) = (PPIIT-PIFFT) / PLAM2 ELSEIF(NC.EQ.2) THEN *** FINT = dlam/lam**2 [(I2 operator I2) - (F operator F)] *** Leading-order terms CALL DFI2DR(2,RH,PIM, ALAMP1,PLAM,F0PI, F2PI) CALL DFI2DR(2,RH,AMETA,ALAMP0,PLAM,F0ET,F2ET) CALL DFI2DR(2,RH,AMETP,ALAMP0,PLAM,F0EP,F2EP) CALL I2NYM(RH,AMR1,ALAMV1,PLAM,F0RO1) CALL I2NYM(RH,AMR2,ALAMV1,PLAM,F0RO2) F0RO=ARO*F0RO1+BRO*F0RO2 CALL I2NYM(RH,A0M,ALAMS1,PLAM,F0A0) CALL OPSEU(RH,F2PI(1),F2PI(2),F2ET(1),F2ET(2),ETFFC,ETFFS,ETFFT) CALL OPSEU(RH,F2PI(1),F2PI(2),F2EP(1),F2EP(2),EPFFC,EPFFS,EPFFT) CALL OSCL(RH,F2PI(1),F2PI(2),F0A0,A0FFS,A0FFT) CALL OVC1(RH,F2PI(1),F2PI(2),F0RO,RFFEES,RFFEET) cp F2A2D=FDEXP(-PLAM2/(4D0*AMA2D2))*I2A2D cp CALL OSCL(RH,F2PI(1),F2PI(2),F2A2D,A2FFS,A2FFT) FINT(1) = (ETIIS-ETFFS) / PLAM2 FINT(2) = (ETIIT-ETFFT) / PLAM2 FINT(3) = (EPIIS-EPFFS) / PLAM2 FINT(4) = (EPIIT-EPFFT) / PLAM2 FINT(5) = (A0IIS-A0FFS) / PLAM2 FINT(6) = (A0IIT-A0FFT) / PLAM2 FINT(7) = (RIIEES-RFFEES) / PLAM2 FINT(8) = (RIIEET-RFFEET) / PLAM2 cp FINT(9) = (A2IIS-A2FFS) / PLAM2 cp FINT(10)= (A2IIT-A2FFT) / PLAM2 ELSEIF(NC.EQ.3) THEN *** FINT = dlam/lam**2 [(I2 operator I2) - (F operator F)] *** 1/M**2 for (1+\kappa) and (1+\kappa)**2 in vector CALL DFI2DR(2,RH,PIM, ALAMP1,PLAM,F0PI, F2PI) CALL DFI2DR(2,RH,AMR1,ALAMV1,PLAM,F0RO1,F2RO1) CALL DFI2DR(2,RH,AMR2,ALAMV1,PLAM,F0RO2,F2RO2) F2RO(1)=ARO*F2RO1(1)+BRO*F2RO2(1) F2RO(2)=ARO*F2RO1(2)+BRO*F2RO2(2) CALL DFI2DR(2,RH,AMOM,ALAMV0,PLAM,F0OM,F2OM) CALL DFI2DR(2,RH,AMFI,ALAMV0,PLAM,F0FI,F2FI) ** Vector-meson functions for 1+f/g and (1+f/g)**2 CALL SPECV(1,RH,PLAM,F2RO,FROK,FROK2,AKRO) CALL SPECV(2,RH,PLAM,F2OM,FOMK,FOMK2,AKOM) CALL SPECV(3,RH,PLAM,F2FI,FFIK,FFIK2,AKFI) ** Operators for 1/M^2 corrections for (1+f/g) and (1+f/g)**2 CALL OVC2(RH,F2PI(1),F2PI(2),FROK(1),FROK2(1),FROK2(2), . RFFEMT,RFFMMS,RFFMMT) CALL OVC3(RH,F2PI(1),F2PI(2),FROK(1),FROK2(1),FROK2(2), . RFFEML,RFFMMC) CALL OVC2(RH,F2PI(1),F2PI(2),FOMK(1),FOMK2(1),FOMK2(2), . OFFEMT,OFFMMS,OFFMMT) CALL OVC2(RH,F2PI(1),F2PI(2),FFIK(1),FFIK2(1),FFIK2(2), . FFFEMT,FFFMMS,FFFMMT) FINT(1) = (RIIEMT-RFFEMT) / PLAM2 FINT(2) = (RIIEML-RFFEML) / PLAM2 FINT(3) = (RIIMMC-RFFMMC) / PLAM2 FINT(4) = (RIIMMS-RFFMMS) / PLAM2 FINT(5) = (RIIMMT-RFFMMT) / PLAM2 FINT(6) = (OIIEMT-OFFEMT) / PLAM2 FINT(7) = (OIIMMS-OFFMMS) / PLAM2 FINT(8) = (OIIMMT-OFFMMT) / PLAM2 FINT(9) = (FIIEMT-FFFEMT) / PLAM2 FINT(10)= (FIIMMS-FFFMMS) / PLAM2 FINT(11)= (FIIMMT-FFFMMT) / PLAM2 ELSEIF(NC.EQ.4) THEN *** FINT = dlam/lam**2 [(I2 operator I2) - (F operator F)] *** Remaining 1/M**2 terms in vector CALL DFI2DR(4,RH,PIM, ALAMP1,PLAM,F0PI,F2PI4) CALL DFI2DR(2,RH,AMR1,ALAMV1,PLAM,F0RO1,F2RO1) CALL DFI2DR(2,RH,AMR2,ALAMV1,PLAM,F0RO2,F2RO2) F0RO=ARO*F0RO1+BRO*F0RO2 F2RO(1)=ARO*F2RO1(1)+BRO*F2RO2(1) F2RO(2)=ARO*F2RO1(2)+BRO*F2RO2(2) CALL I2NYM(RH,AMOM,ALAMV0,PLAM,F0OM) CALL I2NYM(RH,AMFI,ALAMV0,PLAM,F0FI) ** Vector-meson function for 1+2f/g CALL SPECVR(RH,PLAM,F2RO,FRO2K,AKRO) CALL OM2VEC(RH,F2PI4(1),F2PI4(2),FRO2K(1),FRO2K(2), . ROFFM2S,ROFFM2T,ROFFM2L) CALL OM2VK(RH,F2PI4(1),F2PI4(2),F2PI4(3),F2PI4(4), . F0RO,F2RO(1),F2RO(2), ROFFS,ROFFT) CALL OM2VK1(RH,F2PI4(1),F2PI4(2),F2PI4(3),F2PI4(4),F0RO, . ROFFS1,ROFFT1) CALL OM2VK1(RH,F2PI4(1),F2PI4(2),F2PI4(3),F2PI4(4),F0OM, . OMFFS1,OMFFT1) CALL OM2VK1(RH,F2PI4(1),F2PI4(2),F2PI4(3),F2PI4(4),F0FI, . FIFFS1,FIFFT1) FINT(1) = (ROIIM2S-ROFFM2S) / PLAM2 FINT(2) = (ROIIM2T-ROFFM2T) / PLAM2 FINT(3) = (ROIIM2L-ROFFM2L) / PLAM2 FINT(4) = (ROIIS -ROFFS ) / PLAM2 FINT(5) = (ROIIT -ROFFT ) / PLAM2 FINT(6) = (ROIIS1-ROFFS1) / PLAM2 FINT(7) = (ROIIT1-ROFFT1) / PLAM2 FINT(8) = (OMIIS1-OMFFS1) / PLAM2 FINT(9) = (OMIIT1-OMFFT1) / PLAM2 FINT(10)= (FIIIS1-FIFFS1) / PLAM2 FINT(11)= (FIIIT1-FIFFT1) / PLAM2 ELSEIF(NC.EQ.5) THEN *** FINT = dlam/lam**2 [(I2 operator I2) - (F operator F)] *** 1/M**2 terms in scalar and diffractive CALL DFI2DR(3,RH,PIM,ALAMP1,PLAM,F0PI,F2PI3) CALL DFI2DR(2,RH,A0M,ALAMS1,PLAM,F0A0,F2A0) CALL OM2SCA(RH,F2PI3(1),F2PI3(2),F2PI3(3),F2A0(1),F2A0(2), . A0FFM2S,A0FFM2T,A0FFM2L) cp CALL OM2SCA(RH,F2PI3(1),F2PI3(2),F2PI3(3),DI2A2D,DDI2A2D, cp . A2FFM2S,A2FFM2T,A2FFM2L) cp FACA2D=FDEXP(-PLAM2/(4D0*AMA2D2)) FINT(1) = (A0IIM2S-A0FFM2S) / PLAM2 FINT(2) = (A0IIM2T-A0FFM2T) / PLAM2 FINT(3) = (A0IIM2L-A0FFM2L) / PLAM2 cp FINT(4) = (A2IIM2S-A2FFM2S*FACA2D) / PLAM2 cp FINT(5) = (A2IIM2T-A2FFM2T*FACA2D) / PLAM2 cp FINT(6) = (A2IIM2L-A2FFM2L*FACA2D) / PLAM2 ENDIF RETURN END ************************************************************************ SUBROUTINE OPSEU(R,F1,F2,G1,G2,CEN,SIG,TEN) IMPLICIT REAL*8 (A-H,O-Z) C* C* -> (k1.k2)^2 F(r)G(r) = CEN C* -> (sig1.k1xk2)(sig2.k1xk2) F(r)G(r) = C* = SIG (sig1.sig2) + TEN S_{12} C* CEN = F2*G2+2D0*F1*G1/(R*R) SIG = 2D0*(F1*G1/R+F1*G2+F2*G1)/(R*3D0) TEN = ((F1/R-F2)*G1+F1*(G1/R-G2))/(R*3D0) RETURN END ************************************************************************ SUBROUTINE OSCL(R,F1,F2,G0,SIG,TEN) IMPLICIT REAL*8 (A-H,O-Z) C* C* -> (sig1.k1)(sig2.k1) F(r)G(r) = C* = SIG (sig1.sig2) + TEN S_{12} C* SIG =-(F2+2D0*F1/R)*G0/3D0 TEN =-(F2-F1/R)*G0/3D0 RETURN END ************************************************************************ SUBROUTINE OVC1(R,F1,F2,G0,SIG,TEN) IMPLICIT REAL*8 (A-H,O-Z) C* C* -> -(sig1.k1)(sig2.k1) F(r)G(r) = C* = SIG (sig1.sig2) + TEN S_{12} C* SIG = (F2+2D0*F1/R)*G0/3D0 TEN = (F2-F1/R)*G0/3D0 RETURN END ************************************************************************ SUBROUTINE OVC2(R,F1,F2,G1,H1,H2,TEM,SMM,TMM) IMPLICIT REAL*8 (A-H,O-Z) C* C* -> [-k1^2(sig1.k1 sig2.k2 + sig1.k2 sig2.k1) C* +2 k1.k2(sig1.k1)(sig2.k1) ] F(r)G(r) = C* = TEM S_{12} C* -> [ (sig1.k1)(sig2.k1)k2^2 - k1^2 k2^2(sig1.sig2) C* +(sig1.k1xk2)(sig2.k1xk2) ] F(r)G(r) = C* = SMM (sig1.sig2) + TMM S_{12} C* TEM =-2D0*(F2-F1/R)*G1/R SMM =-((F2+F1/R)*(H2+H1/R)+2D0*F1*H1/(R*R))*2D0/3D0 TMM = (2D0*(F2-F1/R)*H2-F2*(H2-H1/R))/3D0 RETURN END ************************************************************************ SUBROUTINE OVC3(R,F1,F2,G1,H1,H2,SOEM,CMM) IMPLICIT REAL*8 (A-H,O-Z) C* C* -> -2i Q.k1xk2(sig1+sig2).k1 F(r)G(r) = SOEM L.S C* -> [ k1^2 k2^2 - (k1.k2)^2 ] F(r)G(r) = CMM C* SOEM= 4D0*F1*G1/(R*R) CMM = 2D0*(F1*H2+F2*H1+F1*H1/R)/R RETURN END ************************************************************************ SUBROUTINE ONAPS(R,F1,F2,F3,G1,G2,G3,CEN,SIG,TEN,SO) IMPLICIT REAL*8 (A-H,O-Z) C* C* -> (k1.k2)^3 F(r)G(r) = CEN C* -> (k1.k2)(sig1.k1xk2)(sig2.k1xk2) F(r)G(r) = C* = SIG (sig1.sig2) + TEN S_{12} C* -> i(k1.k2)(sig1+sig2).(k1xk2) Q.(k1-k2) F(r)G(r) = SO L.S C* R2=R*R CEN =-6D0*(F2-F1/R)*(G2-G1/R)/R2-F3*G3 SIG = 2D0*((F1/R-F2+R*F3)*(G1/R-G2+R*G3)/R2-F3*G3)/3D0 TEN =-((2D0*(F1/R-F2)+0.5D0*R*F3)*(2D0*(G1/R-G2)+0.5D0*R*G3)/R2 . -0.25D0*F3*G3)/3D0 SO =-4D0*(F2-F1/R)*(G2-G1/R)/R2 RETURN END ************************************************************************ SUBROUTINE ONAVS(R,F1,F2,F3,G1,SIG,TEN) IMPLICIT REAL*8 (A-H,O-Z) C* C* -> (k1.k2)(sig1.k1)(sig2.k1) F(r)G(r) = C* = SIG (sig1.sig2) + TEN S_{12} C* SIG = (F3+2D0*(F2-F1/R)/R)*G1/3D0 TEN = (F3-(F2-F1/R)/R)*G1/3D0 RETURN END ************************************************************************ SUBROUTINE OM2VEC(R,F1,F2,G1,G2,SIG,TEN,SO) IMPLICIT REAL*8 (A-H,O-Z) C* C* -> [ -1/2 k1^2(sig1.k1 sig2.k2 + sig1.k2 sig2.k1) C* + (k2^2+k1.k2)(sig1.k1)(sig2.k1) ] F(r)G(r) = C* = SIG (sig1.sig2) + TEN S_{12} C* -> - i Q.k1xk2(sig1+sig2).k1 F(r)G(r) = SO L.S C* SIG = (F2+2D0*F1/R)*(G2+2D0*G1/R) /3D0 TEN = (F2-F1/R)*(G2-G1/R) /3D0 SO = 2D0*F1*G1/(R*R) RETURN END ************************************************************************ SUBROUTINE OM2VK(R,F1,F2,F3,F4,G0,G1,G2,SIG,TEN) IMPLICIT REAL*8 (A-H,O-Z) DATA IPV/1/,IPIMSQ/1/ C* IPIMSQ=1 includes 1/M**2 terms in pseudovector vertices C* C* -> [( k1^2+3k1.k2+k2^2 )(sig1.k1)(sig2.k1) C* - 1/2 k1.k2(sig1.k1 sig2.k2 + sig1.k2 sig2.k1) C* + k1.k2(sig1.k1)(sig2.k1) ] F(r)G(r) = C* = SIG (sig1.sig2) + TEN S_{12} C* SIG = ( (F4+4D0*F3/R)*G0 + (F2+2D0*F1/R)*(G2+2D0*G1/R) . +(3+IPIMSQ)*(F3+2D0*(F2-F1/R)/R)*G1 . - IPIMSQ*IPV*(F2*G2+2D0*F1*G1/(R*R)) )/3D0 TEN = ( (F4+(F3-6D0*(F2-F1/R)/R)/R)*G0 + (F2-F1/R)*(G2+2D0*G1/R) . +(3+IPIMSQ)*(F3-(F2-F1/R)/R)*G1 . - IPIMSQ*IPV*(F2*G2-F1*G1/(R*R)) )/3D0 RETURN END ************************************************************************ SUBROUTINE OM2VK1(R,F1,F2,F3,F4,G0,SIG,TEN) IMPLICIT REAL*8 (A-H,O-Z) C* C* -> k1^2(sig1.k1)(sig2.k1) F(r)G(r) = SIG (sig1.sig2) + TEN S_{12} C* SIG = (F4+4D0*F3/R)*G0/3D0 TEN = (F4+(F3-6D0*(F2-F1/R)/R)/R)*G0/3D0 RETURN END ************************************************************************ SUBROUTINE OM2SCA(R,F1,F2,F3,G1,G2,SIG,TEN,SO) IMPLICIT REAL*8 (A-H,O-Z) DATA IPV/1/,IPIMSQ/1/ C* IPIMSQ=1 includes 1/M**2 terms in pseudovector vertices C* C* -> [ (k2^2+2k1.k2)(sig1.k1)(sig2.k1) C* - 1/2 k1^2(sig1.k1 sig2.k2 + sig1.k2 sig2.k1) C* + 1/2 k1.k2(sig1.k1 sig2.k2 + sig1.k2 sig2.k1) C* - k1.k2(sig1.k1)(sig2.k1) ] F(r)G(r) = C* = SIG (sig1.sig2) + TEN S_{12} C* -> - i Q.k1xk2(sig1+sig2).k1 F(r)G(r) = SO L.S C* SIG = ( (F2+2D0*F1/R)*(G2+2D0*G1/R) . + (1-IPIMSQ)*(F3+2D0*(F2-F1/R)/R)*G1 . + IPIMSQ*IPV*(F2*G2+2D0*F1*G1/(R*R)) ) /3D0 TEN = ( (F2-F1/R)*(G2-G1/R) + (1-IPIMSQ)*(F3-(F2-F1/R)/R)*G1 . + IPIMSQ*IPV*(F2*G2-F1*G1/(R*R)) ) /3D0 SO = 2D0*F1*G1/(R*R) RETURN END ************************************************************************ SUBROUTINE SPECVR(RH,PLAM,D2RO,DRO2K,KAPPA) ** (1+2f/g) function for rho exchange ** Relevant when vector(g) and tensor(f) form factors are different IMPLICIT REAL*8 (A-H,O-Z) REAL*8 KAPPA,D2RO(2),DRO2K(2),FRO1(2),FRO2(2) COMMON/YUKEFF/ARO,AMR1,AROC,AMRC1,AWRC1,BRO,AMR2,BROC,AMRC2,AWRC2, . AVSC1,AVSC2, AEPS,AME1,BEPS,AME2 COMMON/CUTTEN/ALAMT1,ALAVT1,ALAMT0,ALAVT0 CALL DFI2DR(2,RH,AMR1,ALAVT1,PLAM,F0RO1,FRO1) CALL DFI2DR(2,RH,AMR2,ALAVT1,PLAM,F0RO2,FRO2) DRO2K(1)=D2RO(1)+2D0*KAPPA*(ARO*FRO1(1)+BRO*FRO2(1)) DRO2K(2)=D2RO(2)+2D0*KAPPA*(ARO*FRO1(2)+BRO*FRO2(2)) RETURN END ************************************************************************ SUBROUTINE SPECV(N,RH,PLAM,D2VV,D2VT,D2TT,KAPPA) ** (1+f/g) and (1+f/g)**2 functions for rho/omega/phi exchange ** Relevant when vector(g) and tensor(f) form factors are different IMPLICIT REAL*8 (A-H,O-Z) REAL*8 KAPPA,D2VV(2),D2VT(1),D2TT(2), DUM1(2),DUM2(2), . FRO1(2),FRO2(2) COMMON/MESONM/AMPI,AMETA,AMETP, AMRO,AMOM,AMFI, AMA0,AMEP,AMF0, . AMPIC,AMROC,AMA0C,AWPIC,AWROC,AWA0C, AVSC COMMON/YUKEFF/ARO,AMR1,AROC,AMRC1,AWRC1,BRO,AMR2,BROC,AMRC2,AWRC2, . AVSC1,AVSC2, AEPS,AME1,BEPS,AME2 COMMON/CUTTEN/ALAMT1,ALAVT1,ALAMT0,ALAVT0 IF(N.EQ.1) THEN ** rho exchange CALL DFI2DR(2,RH,AMR1,ALAVT1,PLAM,F0RO1,FRO1) CALL DFI2DR(2,RH,AMR2,ALAVT1,PLAM,F0RO2,FRO2) DUM1(1)=ARO*FRO1(1)+BRO*FRO2(1) DUM1(2)=ARO*FRO1(2)+BRO*FRO2(2) CALL DFI2DR(2,RH,AMR1,ALAMT1,PLAM,F0RO1,FRO1) CALL DFI2DR(2,RH,AMR2,ALAMT1,PLAM,F0RO2,FRO2) DUM2(1)=ARO*FRO1(1)+BRO*FRO2(1) DUM2(2)=ARO*FRO1(2)+BRO*FRO2(2) ELSEIF(N.EQ.2) THEN ** omega exchange CALL DFI2DR(2,RH,AMOM,ALAVT0,PLAM,F0OM,DUM1) CALL DFI2DR(2,RH,AMOM,ALAMT0,PLAM,F0OM,DUM2) ELSEIF(N.EQ.3) THEN ** phi exchange CALL DFI2DR(2,RH,AMFI,ALAVT0,PLAM,F0FI,DUM1) CALL DFI2DR(2,RH,AMFI,ALAMT0,PLAM,F0FI,DUM2) ELSE ** STOP '*** SPECV in PIMES: incorrect call' ENDIF D2VT(1)=D2VV(1)+KAPPA*DUM1(1) D2TT(1)=D2VV(1)+KAPPA*(2D0*DUM1(1)+KAPPA*DUM2(1)) D2TT(2)=D2VV(2)+KAPPA*(2D0*DUM1(2)+KAPPA*DUM2(2)) RETURN END *********************************************************************** SUBROUTINE GAUSS(N,X,W,NT) C* This subroutine calculates the weight-factors and argument C* distributions for a Gaussian integration in (-1 , +1). C* N must be 4,8,12,16, or 20. C* otherwise the weight-factors and distributions will be those C* of a two-point Simpson integration. C* NT=1 : N-points Gauss C* NT=2 : N-points Gauss on (-1,0) and N-points Gauss on (0,+1) C----------------------------------------------------------------------- IMPLICIT REAL*8 (A-H,O-Z) DIMENSION WGAUS1(30), WGAUS2(30), WGAUSS(60), X(50), W(50) EQUIVALENCE (WGAUS1(1),WGAUSS(1)), (WGAUS2(1),WGAUSS(31)) DATA WGAUS1 . /.339981043584856D0, .652145154862546D0, . .861136311594053D0, .347854845137454D0, . .183434642495650D0, .362683783378362D0, . .525532409916329D0, .313706645877887D0, . .796666477413627D0, .222381034453374D0, . .960289856497536D0, .101228536290376D0, . .125233408511469D0, .249147045813403D0, . .367831498998180D0, .233492536538355D0, . .587317954286617D0, .203167426723066D0, . .769902674194305D0, .160078328543346D0, . .904117256370475D0, .106939325995318D0, . .981560634246719D0, .047175336386512D0, . .0950125098376374D0, .1894506104550685D0, . .2816035507792589D0, .1826034150449236D0, . .4580167776572274D0, .1691565193950025D0/ DATA WGAUS2 / . .6178762444026437D0, .1495959888165767D0, . .7554044083550030D0, .1246289712555339D0, . .8656312023878317D0, .0951585116824928D0, . .9445750230732326D0, .0622535239386479D0, . .9894009349916499D0, .0271524594117541D0, . .0765265211334973D0, .1527533871307259D0, . .2277858511416451D0, .1491729864726037D0, . .3737060887154196D0, .1420961093183821D0, . .5108670019508271D0, .1316886384491766D0, . .6360536807265150D0, .1181945319615184D0, . .7463319064601508D0, .1019301198172404D0, . .8391169718222188D0, .0832767415767047D0, . .9122344282513259D0, .0626720483341091D0, . .9639719272779138D0, .0406014298003869D0, . .9931285991850949D0, .0176140071391521D0/ IF(MOD(N,4).NE.0) GOTO 2 K=N/4-1 INIT=2*K*(K+1)+1 M=N/2 DO 1 I=1,M J=INIT-2+I*2 X(M+I)=WGAUSS(J) X(M+1-I)=-WGAUSS(J) W(M+I)=WGAUSS(J+1) W(M+1-I)=WGAUSS(J+1) 1 CONTINUE IF(NT.EQ.1) RETURN DO 10 I=1,N X(N+I)=(X(I)+1D0)/2D0 X(I)=(X(I)-1D0)/2D0 W(N+I)=W(I)/2D0 W(I)=W(N+I) 10 CONTINUE RETURN *** Simpson's rule: only odd N 2 H=2D0/(N-1D0) X(1)=-1D0 W(1)=H/3D0 W23=(2D0*H)/3D0 W43=2D0*W23 NN=(N-1)/2 DO 3 I=1,NN J=2*I X(J)=X(J-1)+H W(J)=W43 X(J+1)=X(J)+H W(J+1)=W23 3 CONTINUE X(N)=1D0 W(N)=W(1) RETURN END *********************************************************************** SUBROUTINE AXIAL(R,ISO,SPIN,VSIG,VTEN,VLS,FI0AX,FI1AX,FI2AX) IMPLICIT REAL*8 (A-H,O-Z) DIMENSION FI(7,3) INTEGER SPIN COMMON/COPLNG/ALPV,THPV,PV1, FPI,FETA,FETP,FPI2,FETA2,FETP2, . ALVE,THV ,GV1, GRO,GOM,GFI, GRO2,GOM2,GFI2, . ALVM, FV1, FRO,FOM,FFI, FRO2,FOM2,FFI2, . ALGS,THGS,GS1, GA0,GEP,GF0, GA02,GEP2,GF02, . GFPRO,GFPOM,GFPFI, GFMRO,GFMOM,GFMFI, . FPIC2,GROC2,FROC2,GFPROC,GFMROC,GA0C2 COMMON/CUTOFF/ALAM1,ALAMP1,ALAMV1,ALAMS1, . ALAM0,ALAMP0,ALAMV0,ALAMS0 COMMON/MEANM/ PIM,ROM,A0M COMMON/SCALIN/AMT,AMPV,AMT24 COMMON/AMCOEF/ALF,REDM, AMY,AMY2,AMYI,AMYI2, AMN,AMN2,AMNI,AMNI2, . AMYPN,AMYMN, AMYPN2,AMYMN2, AY2PN2,AY2MN2, . AMYN,AMYNI,AMYNI2, AYPNI2,AYMNI2 COMMON/HBARC/ HBC DATA SR3/1.7320508075688772935D0/, PI/3.14159265358979323846D0/ DATA AX1/0D0/, THAX/0D0/, ALAX/0.4D0/ DATA AMA1/1260D0/, AMF1/1285D0/, AME/1420D0/ RH=R/HBC c FA1= FPI*(AMA1/AMPV) * DSQRT(AMA1**2-ROM**2)/ROM FA1= FPI*(AMA1/AMPV) AX8= FA1*(4D0*ALAX-1D0)/SR3 COST=DCOS(THAX*PI/180D0) SINT=DSIN(THAX*PI/180D0) FF1 =COST*AX8-SINT*AX1 FE =SINT*AX8+COST*AX1 FA12=(4D0*ISO-3D0)*FA1*FA1 FF12=FF1*FF1 FE2 =FE *FE FAC1=FA12/(AMA1*AMA1) FAC2=FF12/(AMF1*AMF1) FAC3=FE2 /(AME *AME ) CALL REGGEFF(RH,AMA1,ALAMP1,FI(1,1)) CALL REGGEFF(RH,AMF1,ALAMP0,FI(1,2)) CALL REGGEFF(RH,AME ,ALAMP0,FI(1,3)) VSIG=-FA12*(FI(1,1)+AMYNI/6D0*FI(2,1)) . -FF12*(FI(1,2)+AMYNI/6D0*FI(2,2)) . -FE2 *(FI(1,3)+AMYNI/6D0*FI(2,3)) . +(FAC1*FI(2,1)+FAC2*FI(2,2)+FAC3*FI(2,3))/3D0 VTEN=+(FA12*FI(6,1)+FF12*FI(6,2)+FE2 *FI(6,3))*AMYNI/4D0 . +(FAC1*FI(6,1)+FAC2*FI(6,2)+FAC3*FI(6,3)) VLS =-(FA12*FI(4,1)+FF12*FI(4,2)+FE2 *FI(4,3))*AMYNI/2D0 FNONLOC=REDM*(4D0*SPIN-3D0)*5D0/6D0*AMYNI FI0AX=-FNONLOC*(FA12*FI(1,1)+FF12*FI(1,2)+FE2*FI(1,3)) FI1AX= FNONLOC*(FA12*FI(4,1)+FF12*FI(4,2)+FE2*FI(4,3))*RH FI2AX=-FNONLOC*(FA12*FI(2,1)+FF12*FI(2,2)+FE2*FI(2,3)) . -FNONLOC*(FA12*FI(4,1)+FF12*FI(4,2)+FE2*FI(4,3))*2D0 RETURN END *********************************************************************** SUBROUTINE PIOME(R,TYPE,EFAC,VPIS,VPIT) C* This subroutine computes the pion potential from the new C* Nijmegen model with the M/E factor. IMPLICIT REAL*8 (A-H,O-Z) REAL*8 PHI(7),PHIC(7) CHARACTER TYPE*2 COMMON/MESONM/AMPI,AMETA,AMETP, AMRO,AMOM,AMFI, AMA0,AMEP,AMF0, . AMPIC,AMROC,AMA0C,AWPIC,AWROC,AWA0C, AVSC COMMON/COPLNG/ALPV,THPV,PV1, FPI,FETA,FETP,FPI2,FETA2,FETP2, . ALVE,THV ,GV1, GRO,GOM,GFI, GRO2,GOM2,GFI2, . ALVM, FV1, FRO,FOM,FFI, FRO2,FOM2,FFI2, . ALGS,THGS,GS1, GA0,GEP,GF0, GA02,GEP2,GF02, . GFPRO,GFPOM,GFPFI, GFMRO,GFMOM,GFMFI, . FPIC2,GROC2,FROC2,GFPROC,GFMROC,GA0C2 COMMON/CUTOFF/ALAM1,ALAMP1,ALAMV1,ALAMS1, . ALAM0,ALAMP0,ALAMV0,ALAMS0 COMMON/AMCOEF/ALF,REDM, AMY,AMY2,AMYI,AMYI2, AMN,AMN2,AMNI,AMNI2, . AMYPN,AMYMN, AMYPN2,AMYMN2, AY2PN2,AY2MN2, . AMYN,AMYNI,AMYNI2, AYPNI2,AYMNI2 COMMON/HBARC/ HBC RH=R/HBC CALL REGGEFF(RH,AMPI,ALAMP1,PHI) FF=FPI2*EFAC VPIS = FF * PHI(2)/3D0 VPIT = FF * PHI(6) IF(TYPE.EQ.'NP' .OR. TYPE.EQ.'PN') THEN CALL REGGEFF(RH,AMPIC,ALAMP1,PHIC) FF=FPIC2/ALF * EFAC VPIS = VPIS + FF*PHIC(2)/3D0 VPIT = VPIT + FF*PHIC(6) ENDIF RETURN END