Multiple elements feature started.

This commit is contained in:
2023-01-24 09:49:19 +01:00
parent 84d97b95c9
commit 12186a8f2f
2 changed files with 37 additions and 51 deletions

View File

@@ -14,7 +14,7 @@ WARN = # -Wall -Wextra
DIALECT = -std=gnu
prefix = /usr/local
# OPS = -c $(DIALECT) $(WARN) $(DEBUG)
FCFLAGS = $(DIALECT) $(WARN) $(DEBUG) -O3 -mcmodel=medium
FCFLAGS = $(DIALECT) $(WARN) $(DEBUG) -O3 -mcmodel=large
all : trimspNL

View File

@@ -50,7 +50,7 @@ c
C These parameters are related to the maximum number of layers MAXNL
C and maximum number of points in the depth distribution MAXD
INTEGER MAXD,MAXD1,MAXD2,MAXD5,MAXDNL5,MAXD2p2,MAXD2meagb
& ,MAXD2meab
& ,MAXD2meab,MAXEL
INTEGER MAXNL,MAXNL5,MAXNLp25,MAXNL5p2,MAXNLm15
C This is the only point where the number of layers and depth
C profile are changed. All other parameters should be changed
@@ -58,18 +58,18 @@ C accordingly.
PARAMETER (MAXD=1000)
PARAMETER (MAXNL=100)
C Maximum number of elements in each layer, was limited to 5.
C PARAMETER (MAXEL=20)
PARAMETER (MAXEL=20)
PARAMETER (MAXD1=MAXD+1)
PARAMETER (MAXD2=MAXD+2)
PARAMETER (MAXD5=MAXD*5)
PARAMETER (MAXD5=MAXD*MAXEL)
PARAMETER (MAXD2p2=MAXD2*MAXD2)
PARAMETER (MAXD2meagb=MAXD2*36*22)
PARAMETER (MAXD2meab=MAXD2*22)
PARAMETER (MAXNL5=MAXNL*5)
PARAMETER (MAXNL5=MAXNL*MAXEL)
PARAMETER (MAXNLp25=MAXNL*MAXNL5)
PARAMETER (MAXDNL5=MAXNL*MAXD5)
PARAMETER (MAXNL5p2=MAXNL5*MAXNL5*MAXD)
PARAMETER (MAXNLm15=(MAXNL-1)*5)
PARAMETER (MAXNLm15=(MAXNL-1)*MAXEL)
LOGICAL TEST(64),TESTR(2000),TEST1(2000)
LOGICAL EQUAL,FORT33
INTEGER*4 ISRCHEQ,ISRCHFGT,ISRCHFGE,ILLZ
@@ -77,7 +77,7 @@ C PARAMETER (MAXEL=20)
INTEGER I,J,IV
INTEGER tryE,negE
INTEGER depth_interval_flag
INTEGER OldNew
INTEGER nel(MAXEL)
INTEGER*4 seconds_start_total,seconds_stop_total
INTEGER*4 NREC1,NREC2,NE1,K,NGIK,ICW
INTEGER*4 ISEED,ISEED2,ISEED3
@@ -115,7 +115,7 @@ C PARAMETER (MAXEL=20)
INTEGER*4 ICSUM,ICSUMS,ICDI,ISPA,ISPAT
INTEGER*4 KK0,KK0R,KK2,KKR,KDEE1,KDEE2
INTEGER*4 NE,NA,NG,NA1,NG1
INTEGER*4 LMAX,JMAX,LJ,INEL,IH,IH1,IY,IY2,IY3
INTEGER*4 JMAX,LJ,INEL,IH,IH1,IY,IY2,IY3
INTEGER*4 JL,KK1,IVMIN,IVMAX,NPA,IREC1,IREC,MAXA,NALL,NSA,KIS
INTEGER*4 IA,IAG,IAGS,IG,IESP,IESLOG
INTEGER*4 IPOT,IPOTR,IRL,ICDIR,ICSBR,ICSUMR,KOI,IGG,KIST
@@ -186,9 +186,9 @@ C REAL Variables
C ZT - atomin numbers, MT - mass numbers (amu), CO - concentration (stoichiometry)
C SBE - surface binding energy, ED - displacement energy, BE - bulk binding energy
C COM - ??
REAL*8 ZT(MAXNL,5),MT(MAXNL,5),CO(MAXNL,5)
& ,SBE(MAXNL,5),ED(MAXNL,5),BE(MAXNL,5),
& COM(5,MAXNL)
REAL*8 ZT(MAXNL,MAXEL),MT(MAXNL,MAXEL),CO(MAXNL,MAXEL)
& ,SBE(MAXNL,MAXEL),ED(MAXNL,MAXEL),BE(MAXNL,MAXEL),
& COM(MAXEL,MAXNL)
REAL*8 MU(MAXNL5,MAXNL5),EC(MAXNL5,MAXNL5),A(MAXNL5,MAXNL5)
& ,F(MAXNL5,MAXNL5) ,KL(MAXNL5,MAXNL5),KOR(MAXNL5,MAXNL5)
& ,KLM(MAXNL,MAXNL5)
@@ -403,47 +403,46 @@ C Ordered as: Number of particles, seed, seed, seed, initial depth, RD, dept
C Third line: Number of layers
read(11,66) NLayers
66 format(9x,I4)
LMAX=NLayers
C Here we read the NLayers structure
DO I=1,NLayers
C Chemical formula (chem), Thickness (DX), density (RHO), and correction factor
C (CK, it is always 1.0??)
C (CK, it is always 1.0??), number of elements (nel)
C Possibly add the number of elements in the layer
READ(11,*) chem(I),DX(I),RHO(I),CK(I)
READ(11,*) chem(I),DX(I),RHO(I),CK(I),nel(I)
C Replace the following lines with read(11,*) (a(i,j), j = 1, 5)
C Atomic numbers
C READ(11,*) ZT(I,1),ZT(I,2),ZT(I,3),ZT(I,4),ZT(I,5)
read(11,*) (ZT(I,j), j = 1, 5)
read(11,*) (ZT(I,j), j = 1, nel(I))
C Mass numbers (amu)
C READ(11,*) MT(I,1),MT(I,2),MT(I,3),MT(I,4),MT(I,5)
read(11,*) (MT(I,j), j = 1, 5)
read(11,*) (MT(I,j), j = 1, nel(I))
C Concentration (stoichiometry)
C READ(11,*) CO(I,1),CO(I,2),CO(I,3),CO(I,4),CO(I,5)
read(11,*) (CO(I,j), j = 1, 5)
read(11,*) (CO(I,j), j = 1, nel(I))
C Surface binding energy
C READ(11,*) SBE(I,1),SBE(I,2),SBE(I,3),SBE(I,4),SBE(I,5)
read(11,*) (SBE(I,j), j = 1, 5)
read(11,*) (SBE(I,j), j = 1, nel(I))
C Displacement energy, always 30 eV??
C READ(11,*) ED(I,1),ED(I,2),ED(I,3),ED(I,4),ED(I,5)
read(11,*) (ED(I,j), j = 1, 5)
read(11,*) (ED(I,j), j = 1, nel(I))
C Bulk binding energy, usually zero
C READ(11,*) BE(I,1),BE(I,2),BE(I,3),BE(I,4),BE(I,5)
read(11,*) (BE(I,j), j = 1, 5)
read(11,*) (BE(I,j), j = 1, nel(I))
C value A-1 of the ziegler tables
C READ(11,*) CH1(I,1),CH1(I,2),CH1(I,3),CH1(I,4),CH1(I,5)
read(11,*) (CH1(I,j), j = 1, 5)
read(11,*) (CH1(I,j), j = 1, nel(I))
C value A-2 of the ziegler tables
C READ(11,*) CH2(I,1),CH2(I,2),CH2(I,3),CH2(I,4),CH2(I,5)
read(11,*) (CH2(I,j), j = 1, 5)
read(11,*) (CH2(I,j), j = 1, nel(I))
C value A-3 of the ziegler tables
C READ(11,*) CH3(I,1),CH3(I,2),CH3(I,3),CH3(I,4),CH3(I,5)
read(11,*) (CH3(I,j), j = 1, 5)
read(11,*) (CH3(I,j), j = 1, nel(I))
C value A-4 of the ziegler tables
C READ(11,*) CH4(I,1),CH4(I,2),CH4(I,3),CH4(I,4),CH4(I,5)
read(11,*) (CH4(I,j), j = 1, 5)
read(11,*) (CH4(I,j), j = 1, nel(I))
C value A-5 of the ziegler tables
C READ(11,*) CH5(I,1),CH5(I,2),CH5(I,3),CH5(I,4),CH5(I,5)
read(11,*) (CH5(I,j), j = 1, 5)
read(11,*) (CH5(I,j), j = 1, nel(I))
ENDDO
1359 CLOSE(UNIT=11)
@@ -497,10 +496,8 @@ C CALCULATION OF CHARGE AND MASS DEPENDENT CONSTANTS
PI2=2.D0*PI
ABC=AB*FP
C This function is used only once to get the number of defined
C layers. It should be replaced!
L=ISRCHEQ(LMAX,DX(1),1,0.D0)-1
L=NLayers
C Checks wether depth interval is an integer denominator of layer thickness or not
C If not, calculated implantation profile is not correct.
C CW is depth increment and DX(K) is the thickness of layer K
@@ -514,16 +511,20 @@ C If the thickness of layer K is not an integer multiples of depth increment
ENDDO
44 CONTINUE
DO I=1,L
C DO I=1,L
C This loop finds out how many elements are in the layer.
C I can be removed in the number is know from the input file.
DO J=1,JMAX
IF(EQUAL(CO(I,J),0.D0)) GOTO 156
ENDDO
J=JMAX+1
C DO J=1,JMAX
C IF(EQUAL(CO(I,J),0.D0)) GOTO 156
C ENDDO
C J=JMAX+1
C I am guessing NJ(I) is the number of elements in layer I
156 NJ(I)=J-1
ENDDO
C 156 NJ(I)=J-1
C ENDDO
C This is the number of elements in layer I
NJ(I) = nel(I)
JT(1) = 0
JT(2) = NJ(1)
LJ= NJ(1)+NJ(2)
@@ -4236,21 +4237,6 @@ C
RETURN
END
INTEGER FUNCTION ISRCHEQ(N,ARRAY,INC,TARGT)
IMPLICIT NONE
INTEGER I,N,J,INC
REAL*8 ARRAY(N),TARGT
J=1
IF(INC.LT.0) J=N*(-INC)
DO I=1,N
IF(ARRAY(J).EQ.TARGT) GO TO 200
J=J+INC
ENDDO
200 ISRCHEQ=I
RETURN
END
SUBROUTINE ENERGGAUSS(ISEED2,Esig,Epar,E0)
IMPLICIT NONE
INTEGER*4 ISEED2