C ******************** Uni1 ************************ C C CALLING PROGRAM FOR E-FIELD AND TUNNEL CURRENT COMPUTATIONS FOR C TIP NEAR UNIFORM SEMICONDUCTOR, FOR PLANAR GEOMETRY. C C VERSION 6.1 - FEB/11 C 6.2 - JUNE/11, RELABEL E, EF IN SIG C C CALLS SEMITIP1 VERSION 6.0 FOR SOLVING POISSON'S EQN (ROUTINES C RHOSURF AND RHOBULK BELOW ARE CALLED BY SEMITIP1 TO SUPPLY CHARGE DENSITIES). C C CALLS SEMIRHO VERSION 6.0 FOR EVALUATING BULK CHARGE DENSITIES. C C CALLS SURFRHO VERSION 6.1 FOR EVALUATING SURFACE CHARGE DENSITIES. C C PARAMETER STATEMENT BELOW, AS WELL AS /CD/ COMMON BLOCK, MUST BE C DUPLICATED IN ROUTINES RHOSURF AND RHOBULK, LOCATED AT BOTTOM OF FILE C PARAMETER(NVDIM=2048,NSDIM=2048,NVDIM1=NVDIM+1, &NVDIM2=2048,NSDIM2=20000,NEDIM=50000) C DIMENSION VAC(2,NVDIM),SEM(2,NSDIM),VSINT(2),S(NSDIM),ITMAX(10), &EP(10),BBIAS(1000) CHARACTER*1 ANS C C /SEMI/ COMMON BLOCK USED BY SEMIRHO ROUTINES C C TK=TEMPERATURE C EGAP=BAND GAP C ED=DONOR BINDING ENERGY C EA=ACCEPTOR BINDING ENERGY C ACB=CONDUCTION BAND EFFECTIVE MASS C AVB=AVERAGED VALENCE BAND EFFECTIVE MASS C CD=CONCENTRATION OF DONORS C CA=CONCENTRATION OF ACCEPTORS C IDEG=INDICATOR OF WHETHER DEGENERACY APPLIES C IINV=INDICATOR OF WHETHER INVERSION IS ALLOWED C COMMON/SEMI/TK,EGAP,ED,EA,ACB,AVB,CD,CA,IDEG,IINV C C /SURF/ COMMON BLOCK USED BY SURFRHO ROUTINES C C ISTK=INDICATOR FOR TEMPERATURE DEPENDENCE OF SURFACE CHARGE C TK1=TEMPERATURE (=TK) C EN0=COMBINED CHARGE NEUTRALITY LEVEL C EN=SEPARATE CHARGE NEUTRALITY LEVELS C DENS=DENSITY FOR SURFACE STATE DISTRIBUTIONS C FWHM=WIDTH OF GAUSSIAN DISTRIBUTIONS C ECENT=CENTROID ENERGY FOR GAUSSIAN DISTRIBUTIONS C COMMON/SURF/ISTK,TK1,EN0,EN(2),DENS(2),FWHM(2),ECENT(2) C C /CD/ COMMON BLOCK USED BY RHOSURF AND RHOBULK, BELOW C C EF=FERMI-LEVEL POSITION C ESTART=STARTING ENERGY IN CHARGE DENSITY TABLES C DELE=ENERGY STEP SIZE IN TABLES C NE=NUMBER OF POINTS IN TABLES C RHOBTAB=TABLE OF BULK TOTAL CHARGE DENSITY C RHOSTAB=TABLE OF SURFACE CHARGE DENSITY C COMMON/CD/EF,ESTART,DELE,NE,RHOBTAB(NEDIM),RHOSTAB(NEDIM) DATA EPSIL0/8.854185E-12/E/1.60210E-19/ PI=4.*ATAN(1.) C C SEMICONDUCTOR AND TIP PARAMETERS C READ(9,*) NPARM DO 900 IPARM=1,NPARM READ(9,*) SEP READ(9,*) CPot READ(9,*) CD READ(9,*) CA WRITE(6,*) ' ' WRITE(16,*) ' ' WRITE(6,*) 'CONTACT POTENTIAL =',CPot WRITE(16,*) 'CONTACT POTENTIAL =',CPot WRITE(6,*) 'DOPING =',CD,CA WRITE(16,*) 'DOPING =',CD,CA READ(9,*) EGAP READ(9,*) ED READ(9,*) EA READ(9,*) ACB read(9,*) AVBH read(9,*) AVBL AVB=exp(2.*alog(sqrt(AVBH**3)+sqrt(AVBL**3))/3.) read(9,*) AVBSO read(9,*) ESO C PARAMETERS BELOW ELIMINATED C read(9,*) E2HH C read(9,*) E2SO READ(9,*) IDEG READ(9,*) IINV IF ((CA.GT.CD.AND.(IINV.EQ.1.OR.IINV.EQ.3)).OR. & (CD.GT.CA.AND.(IINV.EQ.2.OR.IINV.EQ.3))) THEN WRITE(6,*) '****** WARNING - LIKELY INCOMPATIBLE DOPING AND ', & 'INVERSION (IINV) PARAMETER' WRITE(6,*) 'CONTINUE (y/n) ?' READ(5,50) ANS 50 FORMAT(1A1) IF (ANS.NE.'y'.AND.ANS.NE.'Y') STOP END IF READ(9,*) EPSIL READ(9,*) TEM TK=TEM*8.617E-5 TK1=TK READ(9,*) DENS(1) READ(9,*) EN(1) READ(9,*) FWHM(1) READ(9,*) ECENT(1) READ(9,*) DENS(2) READ(9,*) EN(2) READ(9,*) FWHM(2) READ(9,*) ECENT(2) WRITE(6,*) 'FIRST DISTRIBUTION OF SURFACE STATES:' WRITE(16,*) 'FIRST DISTRIBUTION OF SURFACE STATES:' WRITE(6,*) 'SURFACE STATE DENSITY, EN =',DENS(1),EN(1) WRITE(16,*) 'SURFACE STATE DENSITY, EN =',DENS(1),EN(1) WRITE(6,*) 'FWHM, ECENT =',FWHM(1),ECENT(1) WRITE(16,*) 'FWHM, ECENT =',FWHM(1),ECENT(1) WRITE(6,*) 'SECOND DISTRIBUTION OF SURFACE STATES:' WRITE(16,*) 'SECOND DISTRIBUTION OF SURFACE STATES:' WRITE(6,*) 'SURFACE STATE DENSITY, EN =',DENS(2),EN(2) WRITE(16,*) 'SURFACE STATE DENSITY, EN =',DENS(2),EN(2) WRITE(6,*) 'FWHM, ECENT =',FWHM(2),ECENT(2) WRITE(16,*) 'FWHM, ECENT =',FWHM(2),ECENT(2) C C FIND OVERALL CHARGE NEUTRALITY LEVEL, EN0 C IF (DENS(1).EQ.0.AND.DENS(2).EQ.0.) THEN EN0=0. ELSE IF (DENS(1).EQ.0.) THEN EN0=EN(2) ELSE IF (DENS(2).EQ.0.) THEN EN0=EN(1) ELSE CALL ENFIND(EN(1),EN(2),EN0) END IF END IF END IF WRITE(6,*) 'CHARGE-NEUTRALITY LEVEL =',EN0 READ(9,*) ISTK C C GRID SIZES, STEP SIZES, AND ITERATION LIMITS C READ(9,*) NVIN READ(9,*) NSIN READ(9,*) SIZE IF (SIZE.LE.0.) THEN READ(9,*) DELSIN END IF READ(9,*) IPMAX READ(9,*) (ITMAX(IP),IP=1,IPMAX) READ(9,*) (EP(IP),IP=1,IPMAX) READ(9,*) NE C C OUTPUT PARAMETER C READ(9,*) IWRIT C C VOLTAGES C READ(9,*) NBIAS READ(9,*)(BBIAS(IBIAS),IBIAS=1,NBIAS) C C FIND FERMI-LEVEL POSITION C CALL EFFIND(EF) WRITE(6,*) 'FERMI-LEVEL =',EF WRITE(16,*) 'FERMI-LEVEL =',EF RHOCC=RHOCB(EF,0.) RHOVV=RHOVB(EF,0.) WRITE(6,*) 'CARRIER DENSITY IN CB, VB =',RHOCC,RHOVV WRITE(16,*)'CARRIER DENSITY IN CB, VB =',RHOCC,RHOVV C C ****************** LOOP OVER BIAS VOLTAGES **************************** C DO 600 IBIAS=1,NBIAS write(6,*) ' ' write(16,*) ' ' BIAS=BBIAS(IBIAS) write(6,*) 'SEPARATION =',sep write(16,*) 'SEPARATION =',sep C PotTIP=BIAS+CPot write(6,*) ' ' write(16,*) ' ' WRITE(6,*) 'BIAS, TIP POTENTIAL =',BIAS,PotTIP WRITE(16,*) 'BIAS, TIP POTENTIAL =',BIAS,PotTIP C C 1-D SOLUTION FOR BAND BENDING C IF ((CD-CA).EQ.0.) GO TO 105 W=1.E9*SQRT(2.*EPSIL*EPSIL0*AMAX1(1.,ABS(PotTIP))/ & (ABS(CD-CA)*1.E6*E)) WRITE(6,*) '1-D ESTIMATE OF DEPLETION WIDTH (NM) =',W WRITE(16,*) '1-D ESTIMATE OF DEPLETION WIDTH (NM) =',W GO TO 110 105 W=1.E10 C C CONSTRUCT TABLES OF CHARGE DENSITY VALUES C 110 ESTART=AMIN1(EF,EF-PotTIP,EN0) EEND=AMAX1(EF,EF-PotTIP,EN0) ETMP=EEND-ESTART ESTART=ESTART-2.*ETMP EEND=EEND+2.*ETMP DELE=(EEND-ESTART)/FLOAT(NE-1) C PLACE ONE OF THE TABLE VALUES FOR ENERGY AT EF +/- (DELE/2) NETMP=NINT((EF-ESTART)/DELE) ESTART=EF-(NETMP-0.5)*DELE EEND=ESTART+(NE-1)*DELE WRITE(6,*) 'ESTART,EEND,NE =',ESTART,EEND,NE WRITE(16,*) 'ESTART,EEND,NE =',ESTART,EEND,NE WRITE(6,*) 'COMPUTING TABLE OF BULK CHARGE DENSITIES' WRITE(16,*) 'COMPUTING TABLE OF BULK CHARGE DENSITIES' CALL SEMIRHO(DELE,ESTART,NE,NEDIM,RHOBTAB,0,TMP,TMP) WRITE(6,*) 'COMPUTING TABLE OF SURFACE CHARGE DENSITIES' WRITE(16,*) 'COMPUTING TABLE OF SURFACE CHARGE DENSITIES' IF (DENS(1).EQ.0.AND.DENS(2).EQ.0.) THEN DO 115 IE=1,NE RHOSTAB(IE)=0. 115 CONTINUE ELSE CALL SURFRHO(DELE,ESTART,NE,NEDIM,RHOSTAB) END IF C C SEMICONDUCTOR GRID SIZE AND SPACING C NV=NVIN NS=NSIN IF (SIZE.GT.0.) THEN DELS=(W/NS)*SIZE ELSE DELS=DELSIN END IF C C SOLVE THE POTENTIAL PROBLEM C IINIT=1 IWRIT1=IWRIT IF (IWRIT.GT.5) IWRIT1=MOD(IWRIT,5) CALL SEMITIP1(SEP,VAC,SEM,VSINT,S,DELV,DELS,NVDIM,NSDIM, & NV,NS,PotTIP,IWRIT1,ITMAX,EP,IPMAX,Pot0,IERR,IINIT,EPSIL) WRITE(6,120) NS,NV,IERR WRITE(16,120) NS,NV,IERR 120 FORMAT(' RETURN FROM SEMTIP2, NS,NV,IERR =',3I7) WRITE(10,230) SEP,BIAS,CPot,Pot0 230 FORMAT(' ',4G12.4) C C PLOT CROSS-SECTIONAL PROFILE C IF (IWRIT.GE.1) THEN DO 505 J=NV,1,-1 WRITE(11,*) -J*DELV,VAC(1,J) 505 CONTINUE WRITE(11,*) 0.,VSINT(1) DO 510 J=1,NS WRITE(11,*) S(J),SEM(1,J) 510 CONTINUE CLOSE(11) END IF C 550 continue 600 CONTINUE C 900 CONTINUE WRITE(6,*) 'PRESS THE ENTER KEY TO EXIT' WRITE(16,*) 'PRESS THE ENTER KEY TO EXIT' READ(5,*) C stop end C C ROUTINE DEFINING SPECTRUM OF SURFACE STATES C real*8 function sig(ID,ENER) COMMON/SURF/ISTK,TK,EN0,EN(2),DENS(2),FWHM(2),ECENT(2) PI=4.*ATAN(1.) c c sig=0. c if (ENER.lt.0.) return c if (ENER.gt.egap) return c if (fwhm(ID).eq.0.) go to 200 width=fwhm(ID)/(2.*sqrt(2.*alog(2.))) sig=-dexp(-1.d0*(ENER-(EN(ID)+ecent(ID)))**2/(2.*width**2))+ & dexp(-1.d0*(ENER-(EN(ID)-ecent(ID)))**2/(2.*width**2)) sig=sig*dens(ID)/(SQRT(2.*pi)*width) return c c sig=exp(-(ENER-EN)**2/(2.*width**2)) c sig=sig*dens/(SQRT(2.*pi)*width) c c sig=exp(-ENER/width)- c & exp((ENER-egap)/width) c sig=sig*dens/WIDTH c 200 sig=dens(ID) if (ENER.gt.en(ID)) sig=-sig return end C C ROUTINE DEFINING SPATIAL DISTRIBUTION OF SURFACE STATES C (IN THIS 1D CASE, NO SPATIAL VARIATION IS POSSIBLE) C REAL FUNCTION RHOSURF(Pot) C PARAMETER(NVDIM=2048,NSDIM=2048,NVDIM1=NVDIM+1, &NVDIM2=2048,NSDIM2=20000,NEDIM=50000) COMMON/CD/EF,ESTART,DELE,NE,RHOBTAB(NEDIM),RHOSTAB(NEDIM) PI=4.*ATAN(1.) C ENER=EF-Pot IENER=NINT((ENER-ESTART)/DELE)+1 RHO=0. IF (IENER.GE.1.AND.IENER.LE.NE) THEN RHO=RHOSTAB(IENER) ELSE RHO=RHOS(ENER,DELE) END IF C RHOSURF=RHO RETURN END C C ROUTINE DEFINING SPATIAL DISTRIBUTION OF BULK STATES C REAL FUNCTION RHOBULK(Pot,S,J,NS) C PARAMETER(NVDIM=2048,NSDIM=2048,NVDIM1=NVDIM+1, &NVDIM2=2048,NSDIM2=20000,NEDIM=50000) COMMON/CD/EF,ESTART,DELE,NE,RHOBTAB(NEDIM),RHOSTAB(NEDIM) C ENER=EF-Pot IENER=NINT((ENER-ESTART)/DELE)+1 RHO=0. IF (IENER.GE.1.AND.IENER.LE.NE) THEN RHO=RHOBTAB(IENER) ELSE RHO=RHOB(ENER,0.) END IF C C ADD STATEMENT BELOW TO RESTRICT REGION OF BULK CHARGE C IF (S.LT.10.) RHO=0. C RHOBULK=RHO RETURN END