C N. R. BADNELL PROJECT V1.4 20/05/00 C C PROJECT PSEUDO-STATES ONTO PHYSICAL STATES C PROGRAM MAIN CSUB SUBROUTINE PROJECT(EPMIN,EPMAX,SUMN3) C IMPLICIT DOUBLE PRECISION (A-H,O-Z) C PARAMETER(MZLR1=6) !MAX ORBITAL L PARAMETER(MZNR1=50) !MAX PHYSICAL N PARAMETER(MZNR2=25) !MAX PSEUDO N PARAMETER(NTMS=100) PARAMETER(MZTAR1=NTMS) !NO. OF (CC) TARGET TERMS PARAMETER(MZTAR2=300) !NO. OF (CI) TARGET TERMS C PARAMETER (MXORB=((MZLR1+1)*(2*MZNR2-MZLR1))/2) !NO. OF PSEUDO NL C PARAMETER(ZERO=0.0D0) PARAMETER(HALF=0.5D0) PARAMETER(ONE=1.0D0) C DIMENSION OVRLAP(0:MZLR1,MZNR1,MZNR2) DIMENSION INDX2(MZNR2,0:MZLR1) DIMENSION NCONF(MZTAR1),CMIX(MZTAR2,MZTAR1) DIMENSION IQTERM(MXORB,MZTAR2,MZTAR1) DIMENSION NI(MXORB),LI(MXORB) DIMENSION ISTERM(MZTAR1),LTERM(MZTAR1),ETERM(MZTAR1) X,LORDER(MZTAR1),NORDER(MZTAR1) DIMENSION SUMM(MZTAR1,MZNR1),SUMN(MZTAR1),SUMN3(MZTAR1) C NAMELIST/SPROJ/IPRINT,EPMIN,EPMAX !COMMENT-OUT FOR SUBROUTINE C OPEN(5,FILE='pin',STATUS='UNKNOWN') !COMMENT-OUT FOR SUBROUTINE OPEN(6,FILE='pout',STATUS='UNKNOWN') !COMMENT-OUT FOR SUBROUTINE C OPEN(18,FILE="MIXDW.DAT",STATUS="OLD") ! TARGET FROM XTRA "DW" STG2 RUN OPEN(19,FILE='OVRLAP',STATUS='OLD') ! OVERLAPS FROM XTRA AS RUN C IPRINT=-4 C C ONLY LOOK AT PROJECTION ONTO PHYSICAL BOUND FOR PSEUDO STATES LYING C BETWEEN EPMIN AND EPMAX (RYDBERGS RELATIVE TO GROUND STATE) C EPMIN=-1.0 !COMMENT-OUT FOR SUBROUTINE EPMAX=99999.0 !COMMENT-OUT FOR SUBROUTINE C READ(5,SPROJ,END=22) !COMMENT-OUT FOR SUBROUTINE C WRITE(6,311)EPMIN,EPMAX C IF(EPMAX.LT.EPMIN)WRITE(6,*)'WARNING: NO PROJECTION' C C C----BEGIN INPUT OF TARGET E-VECTORS/ENERGIES USED IN R-MATRIX CALCN. C C NTERM--NUMBER OF LS TERMS IN TARGET C MAXORB--NUMBER OF ORBITALS THAT ARE USED TO DEFINE A CONFIGURATION. C KCOR--CLOSED SHELLS DEFINED IN STG2. C 22 READ(18,7,END=400) NTERM,MAXORB,KCOR 7 FORMAT(3I5) C IF(NTERM.GT.MZTAR1)THEN WRITE(6,*)' **ERROR** NTERM = ',NTERM,' > ',MZTAR1 STOP 600 ENDIF IF(MAXORB.GT.MXORB)THEN WRITE(6,*)' **ERROR** MAXORB = ',MAXORB,' > MXORB = ',MXORB STOP 601 ENDIF IF(MAXORB.LE.0)THEN WRITE(6,*)' **ERROR** MAXORB = ',MAXORB,' .LE. 0' STOP 602 ENDIF C C READ THE NL'S DEFINED BY R-MATRIX STG2 C READ(18,*)(NI(J),LI(J),J=1,MAXORB) C C C NTERM--NUMBER OF LS TERMS IN TARGET C ISTERM--PARITY AND SPIN ANGULAR MOMENTUM OF LS TERM C LTERM--ORBITAL ANGULAR MOMENTUM OF LS TERM C ETERM--ENERGY IN RYDBERGS OF LS TERM C NCONF--NUMBER OF CONFIGURATIONS FOR LS TERM C CMIX--MIXING COEFFICIENT FOR EACH CONFIGURATION IN LS TERM C IQTERM--OCCUPATION NUMBER FOR EACH SUBSHELL C FOR EACH CONFIGURATION IN LS TERM C IF(IPRINT.GE.0)WRITE(6,170) NTERM 170 FORMAT(///5X, 'NTERM=',I5/) C DO 71 J=1,NTERM C READ(18,172) ISTERM(J),LTERM(J),NCONF(J),ETERM(J),METAJ 172 FORMAT(3I5,F12.6,I5) C IF(IPRINT.GE.0) XWRITE(6,173) J,ISTERM(J),LTERM(J),NCONF(J),ETERM(J),METAJ 173 FORMAT(5X,I3,3X, '(2S+1)=',I3,3X, 'L=',I3,3X, 'NCONF=',I3,3X, X 'ETERM=',F10.6, ' RY',3X,'META=',I3) C NNCONF=NCONF(J) IF(NNCONF.GT.MZTAR2)THEN WRITE(6,*)' **ERROR** NNCONF = ',NNCONF,' > MZTAR2 = ',MZTAR2 STOP 603 ENDIF C DO 84 I=1,NNCONF C READ(18,*) CMIX(I,J),(IQTERM(II,I,J),II=1,MAXORB) C IF(IPRINT.GT.0) XWRITE(6,186) I,CMIX(I,J1),(IQTERM(II,I,J1),II=1,MAXORB) 186 FORMAT(10X,I3,5X,'CMIX=',F10.6,5X, 'IQTERM=',(20I3)) C 84 CONTINUE 71 CONTINUE C C SET-UP ENERGY ORDER INDEXING FOR ALL OF THE TERMS C N.B. FIRST TERM SHOULD BE THE GROUND (ETERM(1)=0.0) C CALL ORDER(ETERM,NORDER,NTERM,1) C C SET-UP REVERSE INDEX - NOT USED YET. C DO 87 I=1,NTERM LORDER(NORDER(I))=I 87 CONTINUE C C--END E-VECTOR/ENERGY INPUT C C C ZEROISE C DO 1 L=0,MZLR1 DO 2 N=1,MZNR2 DO 3 M=1,MZNR1 OVRLAP(L,M,N)=ZERO 3 CONTINUE 2 CONTINUE 1 CONTINUE C C N1 IS PHYSICAL N, N2 IS PSEUDO PLUS FEW LOW-LYING PHYSICAL. C N1MAX=0 N1MIN=9999 C READ(19,*,END=401)IDUM REWIND(10) C 4 READ(19,*,END=5)LA,N1,N2,DS C IF(N1.GT.MZNR1)THEN WRITE(6,*)' **ERROR**: N1 =',N1,' > MZNR1 = ',MZNR1 STOP 604 ENDIF C N0=IABS(N2) IF(N0.GT.MZNR2)THEN WRITE(6,*)' **ERROR**: N2 =',N0,' > MZNR2 = ',MZNR2 STOP 605 ENDIF C IF(LA.GT.MZLR1)THEN WRITE(6,*)' **ERROR**: LA =',LA,' > MZLR1 = ',MZLR1 STOP 606 ENDIF C INDX2(N0,LA)=N2 N1MAX=MAX(N1MAX,N1) N1MIN=MIN(N1MIN,N1) C OVRLAP(LA,N1,N0)=DS GO TO 4 C 5 CONTINUE C C EVALUATE PROJECTION OF PSEUDO ON TO PHYSICAL DISCRETE C WRITE(6,308) IF(IPRINT.LT.-2.AND.IPRINT.GT.-5)WRITE(6,310) C C--BEGIN LOOP OVER TERMS C DO 50 J=1,NTERM !J IS ENERGY ORDERED C J0=NORDER(J) !J0 IS ORIGINAL ORDER DO 55 N=1,MZNR1 SUMM(J,N)=ZERO 55 CONTINUE C C ONLY LOOK AT STATES LYING BETWEEN EPMIN AND EPMAX: C SUMN3(J)=ZERO !CONTINUUM FRACTION HERE IF(ETERM(J0).LT.EPMIN)GO TO 50 SUMN3(J)=ONE !CONTINUUM FRACTION HERE IF(ETERM(J0).GT.EPMAX)GO TO 50 C NNCONF=NCONF(J0) DO 60 I=1,NNCONF C C SKIP IF MIXING NEGLIGIBLE IF(ABS(CMIX(I,J0)).LT.1.0D-6)GO TO 60 C C LOOK FOR PSEUDO ORBITAL IN CONFIG LIST DO 65 M=1,MAXORB NL=IQTERM(M,I,J0) IF(NL.EQ.0)GO TO 65 !UNOCCUPIED L=LI(M) NB=NI(M) N0=INDX2(NB,L) IF(N0.GT.0)GO TO 65 GO TO 66 !ORBITAL IS PSEUDO 65 CONTINUE C PSEUDO IS MIXED WITH PHYSICAL SO UNIT OVERLAP WITH LAST N SUMM(J,NB)=SUMM(J,NB)+CMIX(I,J0) GO TO 60 C 66 DO 70 N=N1MIN,N1MAX SUMM(J,N)=SUMM(J,N)+OVRLAP(L,N,NB)*CMIX(I,J0) 70 CONTINUE C 60 CONTINUE !END LOOP OVER MIXING C C SQUARE-UP AND SUM OVER ALL PHYSICAL IF(IPRINT.GE.-2)WRITE(6,300)J,ISTERM(J0),LTERM(J0),ETERM(J0) C SUMN(J)=ZERO DO 75 N=1,N1MAX SUMM(J,N)=SUMM(J,N)*SUMM(J,N) SUMN(J)=SUMN(J)+SUMM(J,N) SUMN3(J)=SUMN(J)+HALF*(N-1)*SUMM(J,N) IF(IPRINT.GE.-2)WRITE(6,301)N,SUMM(J,N),SUMN(J),SUMN3(J) 75 CONTINUE C IF(IPRINT.GE.-2)WRITE(6,302)SUMN(J),SUMN3(J) IF(IPRINT.LT.-2.AND.IPRINT.GT.-5) XWRITE(6,309)J,ISTERM(J0),LTERM(J0),ETERM(J0),SUMN(J),SUMN3(J) C C CONVERT TO CONTINUUM FRACTION C SUMN3(J)=ONE-SUMN3(J) C 50 CONTINUE !END LOOP OVER TERMS C C C DO REVERSE PRINTOUT I.E. PHYSICAL ONTO PSEUDO'S C ONLY MEANINGFUL FOR TARGETS ALL WITH SAME ORBITAL L. C IF(IPRINT.GE.-3)THEN WRITE(6,307) DO 80 N=1,N1MAX SUMP=ZERO IF(IPRINT.GE.-1)WRITE(6,303)N DO 90 J=1,NTERM IF(IPRINT.GE.-1)WRITE(6,304)J,SUMM(J,N) SUMP=SUMP+SUMM(J,N) 90 CONTINUE IF(IPRINT.GE.-1)WRITE(6,305)SUMP IF(IPRINT.LT.-1)WRITE(6,306)N,SUMP 80 CONTINUE ENDIF C C STOP !COMMENT-OUT FOR SUBROUTINE CSUB RETURN C 400 WRITE(6,*)'**ERROR** IN PROJECT: MIXDW.DAT FILE IS EMPTY' STOP 607 401 WRITE(6,*)'**ERROR** IN PROJECT: OVRLAP FILE IS EMPTY' STOP 608 C 300 FORMAT(//'TERM=',I3,5X,'(2S+1) L= ',2I2,5X,'ETERM=',F10.6) 301 FORMAT(10X,'N=',I3,5X,'SUMM=',F8.4,5X,'SUMN=',F8.4,5X X,'SUM =',F8.4) 302 FORMAT(/20X,'SUMN=',F8.4,23X,'SUM =',F8.4) 303 FORMAT(//'N=',I3) 304 FORMAT(5X,'T=',I3,5X,'SUMM=',F8.4) 305 FORMAT(/15X,'SUMP=',F8.4) 306 FORMAT('N=',I3,15X,'SUMP=',F8.4) 307 FORMAT(//'PHYSICAL ONTO PSEUDO') 308 FORMAT(//'PSEUDO ONTO PHYSICAL') 309 FORMAT(1X,I3,10X,2I2,11X,F10.6,5X,F8.4,5X,F8.4) 310 FORMAT(//'TERM',6X,'(2S+1) L',15X,'ETERM',9X,'SUMN',9X,'SUM') 311 FORMAT(//'PSEUDO STATES LYING BETWEEN EPMIN AND EPMAX ARE' X,' PROJECTED ONTO THE PHYSICAL CONTINUUM:',2F12.6) C END C C ********************************************************************** C SUBROUTINE ORDER(EN,NORDER,NDIM,IUP) IMPLICIT DOUBLE PRECISION (A-H,O-Z) C C RETURNS NORDER(I)=POINTER TO I-TH ENERGY IN EN ARRAY, C IUP=1 FOR ASCENDING ENERGIES, IUP=-1 FOR DESCENDING ENERGIES C DIMENSION EN(NDIM),NORDER(NDIM) C DO 40 K = 1,NDIM J = K J1 = J - 1 IF (J1.EQ.0) GOTO 30 E = EN(J) + IUP*1.0D-7 DO 20 I = 1,J1 IF (J.LT.K) GOTO 10 N = NORDER(I) IF (IUP.GT.0 .AND. E.GT.EN(N)) GOTO 20 IF (IUP.LT.0 .AND. E.LT.EN(N)) GOTO 20 10 CONTINUE NORDER(J) = NORDER(J-1) J = J - 1 20 CONTINUE 30 NORDER(J) = K 40 CONTINUE C END C C **********************************************************************