C C See eq. 9.15 of Kris Heyde, "The Nuclear Shell Model", Springer-Verlag, C 1994. C PROGRAM DELT REAL NU,NORM1,NORM2 REAL INT1,INT2 INTEGER PAR COMMON/SHELL/E(20),N(20),L(20),J(20),IT(20),NUM(20) COMMON/PARAM/ALFA,VEFF,IDENT COMMON/VN06FC/FCT(40) COMMON/VN07CB/NU,AKNL(4,4,8),ANL(4,8),INT1(20),INT2(20) C C C MATRIX ELEMENTS IN ORDER C (L+S)J IN THAT ORDER C R(N,L) POSITIVE AT THE ORIGIN C N=1,2,3,.... C JA=2*J DOUBLE PHYSICAL VALUE C LA SINGLE PHYSICAL VALUE C V=-VEFF*(1-ALPHA+ALPHA*SIGMA(P)*SIGMA(N))*DELTA FUNCTION C NU=(M*OMEGA)/(2.*HBAR) C C GOOD ESTIMATES ARE VEFF=400;ALPHA=0.2 C IDENT=1 IDENTICAL PARTICLES C C INPUT DATA FOR PN-MATRIX ELEMENTS NU,ALFA,VEFF,IDENT C CALL VN07AL WRITE(3,360)NU 360 FORMAT(1H ,'NU=',F10.8/) READ(1,400)ALFA,VEFF,IDENT 400 FORMAT(2F8.4,I2) WRITE(3,365)ALFA,VEFF 365 FORMAT(1H ,'ALFA=',F8.4,'VEFF=',F8.4) FCT(1)=1. DO 50 I=2,40 FCT(I)=FCT(I-1)*(I-1.)/10. 50 CONTINUE C C IMAX TOTAL NUMBER OF S.P .ORBITS C READ(1,31)IMAX 31 FORMAT(I2) WRITE(3,26)IMAX 26 FORMAT(I4) DO 38 I=1,IMAX READ(1,32)E(I),N(I),L(I),J(I),IT(I),NUM(I) 32 FORMAT(F5.3,5I4) WRITE(3,33)E(I),N(I),L(I),J(I),IT(I),NUM(I) 33 FORMAT(1H ,'E=',F5.3,' N L J=',3I4,'/2 IT=',I2,'NUM=',I4) 38 CONTINUE 60 READ(1,34)I1,I2,I3,I4,JP,ICON 34 FORMAT(6I2) X=DELTA(I1,I2,I3,I4,JP) WRITE(3,35)X 35 FORMAT(1H ,'MATRIX EL.=',E20.8) IF (ICON.EQ.0) GOTO 60 STOP END FUNCTION DELTA(I1,I2,I3,I4,JP) COMMON/PARAM/ALFA,VEFF,IDENT COMMON/SHELL/E(20),N(20),L(20),J(20),IT(20),NUM(20) REAL NORM1,NORM2 NA=N(I1)-1 NB=N(I2)-1 NC=N(I3)-1 ND=N(I4)-1 A1=0.5 IF(IDENT.EQ.1)A1=1. IF(JP.EQ.0.AND.IDENT.NE.1)GOTO 106 IF(JP.GE.1.AND.IDENT.NE.1)GOTO 101 106 ZL=0. GOTO 102 101 ZL=S(J(I2),J(I1),JP)*S(J(I4),J(I3),JP)/(4.*JP*(JP+1.)) 102 XL=SQRT((J(I1)+1.)*(J(I2)+1.)*(J(I3)+1.)*(J(I4)+1.))*(-1)**((J(I1) 1+J(I3))/2+L(I1)+L(I3))/(2.*JP+1.)* 2VN02BA(J(I2),J(I1),2*JP,1,-1,0)* 3VN02BA(J(I4),J(I3),2*JP,1,-1,0)* 4(A1*(1.-(-1)**((J(I2)+J(I4))/2+L(I1)+L(I3))*ZL)- 5ALFA*(1.+(-1)**(L(I1)+L(I2)+JP)))*VEFF* 6SLATER(NA,L(I1),NB,L(I2),NC,L(I3),ND,L(I4)) NORM1=1. NORM2=1. IF(IDENT.NE.1)GOTO 104 IF(I1.EQ.I2)NORM1=1./SQRT(2.) IF(I3.EQ.I4)NORM2=1./SQRT(2.) 104 XL=XL*NORM1*NORM2 DELTA=XL RETURN END FUNCTION S(IA,IB,IC) C HULPFUNKTIE S=IA+1+(IB+1)*((-1)**((IA+IB)/2+IC)) RETURN END FUNCTION VN02BA(J1,J2,J,M1,M2,M) COMMON/VN06FC/FCT(40) INTEGER Z,ZMIN,ZMAX,FASE C VN02BA DENOTES THE CLEBSCH-GORDAN COEFFICIENT C C CC=0. IF(M1+M2-M)20,1,20 1 IF(IABS(M1)-IABS(J1))2,2,20 2 IF(IABS(M2)-IABS(J2))3,3,20 3 IF(IABS(M)-IABS(J))4,4,20 4 IF(J-J1-J2)5,5,20 5 IF(J-IABS(J1-J2))20,6,6 6 ZMIN=0 IF(J-J2+M1)7,8,8 7 ZMIN=-J+J2-M1 8 IF(J-J1-M2+ZMIN)9,10,10 9 ZMIN=-J+J1+M2 10 ZMAX=J1+J2-J IF(J2+M2-ZMAX)11,12,12 11 ZMAX=J2+M2 12 IF(J1-M1-ZMAX)13,14,14 13 ZMAX=J1-M1 14 DO 15 Z=ZMIN,ZMAX,2 JA=Z/2+1 JB=(J1+J2-J-Z)/2+1 JC=(J1-M1-Z)/2+1 JD=(J2+M2-Z)/2+1 JE=(J-J2+M1+Z)/2+1 JF=(J-J1-M2+Z)/2+1 FASE=((-1)**(Z/2)) F2=FASE 15 CC=CC+F2/(FCT(JA)*FCT(JB)*FCT(JC)*FCT(JD)*FCT(JE)* 1FCT(JF)) JA=(J1+J2-J)/2+1 JB=(J1-J2+J)/2+1 JC=(-J1+J2+J)/2+1 JD=(J1+M1)/2+1 JE=(J1-M1)/2+1 JF=(J2+M2)/2+1 JG=(J2-M2)/2+1 JH=(J+M)/2+1 JI=(J-M)/2+1 JJ=(J1+J2+J+2)/2+1 F1=J+1 CC=SQRT(F1*FCT(JA)*FCT(JB)*FCT(JC)*FCT(JD)*FCT(JE)* 1FCT(JF)*FCT(JG)*FCT(JH)*FCT(JI)/FCT(JJ))*CC 20 VN02BA=CC /SQRT(10.) RETURN END FUNCTION SLATER(NI,LI,NII,LII,NK,LK,NKK,LKK) C C FUNCTION SLATER CALCULATES SLATER-INTEGRALS WITH A DELTA-FORCE C USING FOUR RADIAL HARMONIC OSCILLATOR WAVEFUNCTIONS WITH QUANTUM C NUMBERS (NI,LI),(NII,LII),(NK,LK),(NKK,LKK) AND N=0,1,2,3. C C (NI,LI),(NII,LII),(NK,LK),(NKK,LKK): THE FOUR SETS OF QUANTUM C NUMBERS OF THE 4 RADIAL WAVEFUNCTIONS WITH N=0,1,2,3 AND C L=0,1,2,3,4,5,6,7 C INT: THE SLATER-INTEGRAL (DIVISION BY 4*PI INCLUDED) C COMMON/VN07CB/NU,AKNL(4,4,8),ANL(4,8),INT1(20),INT2(20) REAL NU REAL INT1,INT2 ROM=0. L=LI+LII+LK+LKK NI1=NI+1 NII1=NII+1 NK1=NK+1 NKK1=NKK+1 LI1=LI+1 LII1=LII+1 LK1=LK+1 LKK1=LKK+1 L1=(L+1)/2 IF(L/2.EQ.L/2.) GO TO 1 DO 2 I=1,NI1 DO 2 II=1,NII1 DO 2 K=1,NK1 DO 2 KK=1,NKK1 L2=L1+I+II+K+KK-4 2 ROM=ROM+AKNL(I,NI1,LI1)*AKNL(II,NII1,LII1)*AKNL(K,NK1,LK1)*AKNL(KK 1,NKK1,LKK1)*INT1(L2) GO TO 3 1 DO 4 I=1,NI1 DO 4 II=1,NII1 DO 4 K=1,NK1 DO 4 KK=1,NKK1 L2=L/2+I+II+K+KK-3 4 ROM=ROM+AKNL(I,NI1,LI1)*AKNL(II,NII1,LII1)*AKNL(K,NK1,LK1)*AKNL(KK 1,NKK1,LKK1)*INT2(L2) 3 SOM=ROM SLATER=ANL(NI1,LI1)*ANL(NII1,LII1)*ANL(NK1,LK1)*ANL(NKK1,LKK1)*SOM 1*(2*NU)**1.5/78.95683523 RETURN END SUBROUTINE VN07AL C C SUBROUTINE CALCULATES A NUMBER OF COEFFICIENTS NECESSARY FOR THE C CALCULATIONS OF THE SLATER-INTEGRALS C C FAK1(N),N=1,4: FAK1(N)=(N-1)! C FAK2(N),N=1,11: FAK2(N)=(2N-1)!! C ANL(N,L),N=1,4, L=1,8: C ANL(N,L)=SQRT(2**(L-N+2)*(2*NU)**(L+1/2)*(2L+2N-3)!!)/SQRT(SQRT(PI) C *((2L-1)!!)**2*(N-1)!) C AKNL(K,N,L),K=1,4, N=1,4, L=1,8: C AKNL(K,N,L)=(-2)**(K-1)*(N-1)!*(2L-1)!!/(K-1)!*(N-K)!*(2L+2K-3)!!) C INT1(I),I=1,20: INT1(I)=I!/2**(I+1) C INT2(I),I=1,20: INT2(I)=SQRT(PI)*(2I-1)!!/2**(2I+1/2) C NU=M*OMEGA/2*H-BAAR C REAL NU REAL INT1,INT2 COMMON/VN07CB/NU,AKNL(4,4,8),ANL(4,8),INT1(20),INT2(20) DIMENSION FAK1(21),FAK2(20) READ(1,1) NU 1 FORMAT(F10.8) FAK1(1)=1. DO 2 IR=1,20 2 FAK1(IR+1)=FAK1(IR)*IR FAK2(1)=1. DO 3 IR=1,19 3 FAK2(IR+1)=FAK2(IR)*(2*IR+1) DO 4 N=1,4 DO 4 L=1,8 LNS=L+N F1=FAK2(LNS-1)/(FAK2(L)**2*FAK1(N)) ANL(N,L)=SQRT(2.**(L-N+2.)*F1) DO 4 K=1,N KNS=N-K LKS=L+K 4 AKNL(K,N,L)=(-2)**(K-1)*FAK1(N)*FAK2(L)/(FAK1(K)*FAK1(KNS+1)*FAK2( 1LKS-1)) DO 5 I=1,20 INT1(I)=FAK1(I+1)/2**(I+1) 5 INT2(I)=1.25331413731550*FAK2(I)/4.**I RETURN END C C COMMAND FILE FOR THE CALCULATION OF A NUMBER OF MATTRIX ELEMENTS C $ ASSIGN SYS$INPUT FOR001 $ ASSIGN SYS$OUTPUT FOR003 R DELTA.EXE 0.083 0.5000300.0000 2 0.000 1 5 9 1 1 0.000 1 6 13 0 2 1 2 1 2 2 0 1 2 1 2 3 0 1 2 1 2 4 0 1 2 1 2 5 0 1 2 1 2 6 0 1 2 1 2 7 0 1 2 1 2 8 0 1 2 1 2 9 0 1 2 1 210 0 1 2 1 211 1