Commit e60e6434 authored by sjplimp's avatar sjplimp
Browse files

git-svn-id: svn://svn.icms.temple.edu/lammps-ro/trunk@4666 f3b2605a-c512-4ea7-a41b-209d697bcdaa
parent bce3ded4
Loading
Loading
Loading
Loading
+9 −0
Original line number Diff line number Diff line
 &funccard
 atomtype='Ta'
 &end
 &funccard
 atomtype='Cu'
 &end
 &funccard
 &end
+448 −0
Original line number Diff line number Diff line
Cu
2.556162
1.554485
21.175871
21.175395
8.127620
4.334731
0.396620
0.548085
0.308782
0.756515
-2.170269
-0.263788
1.088878
-0.817603
-2.19
0.00
0.561830
-2.100595
0.310490
-2.186568
29
63.546
-2.100595
4.334731
0.756515
0.85
1.15
Ag
2.891814
1.106232
14.604100
14.604144
9.132010
4.870405
0.277758
0.419611
0.339710
0.750758
-1.729364
-0.255882
0.912050
-0.561432
-1.75
0.00
0.744561
-1.150650
0.783924
-1.748423
47
107.8682
-1.150650
4.870405
0.750758
0.85
1.15
Au
2.885034
1.529021
19.991632
19.991509
9.516052
5.075228
0.229762
0.356666
0.356570
0.748798
-2.937772
-0.500288
1.601954
-0.835530
-2.98
0.00
1.706587
-1.134778
1.021095
-2.978815
79
196.96654
-1.134778
5.075228
0.748798
0.85
1.15
Ni
2.488746
2.007018
27.562015
27.562031
8.383453
4.471175
0.429046
0.633531
0.443599
0.820658
-2.693513
-0.076445
0.241442
-2.375626
-2.70
0.00
0.265390
-0.152856
0.445470
-2.7
28
58.6934
-0.152856
4.471175
0.820658
0.85
1.15
Pd
2.750897
1.595417
21.335246
21.940073
8.697397
4.638612
0.406763
0.598880
0.397263
0.754799
-2.321006
-0.473983
1.615343
-0.231681
-2.36
0.00
1.481742
-1.675615
1.130000
-2.352753
46
106.42
-1.675615
4.638612
0.754799
0.85
1.15
Pt
2.771916
2.336509
33.367564
35.205357
7.105782
3.789750
0.556398
0.696037
0.385255
0.770510
-1.455568
-2.149952
0.528491
1.222875
-4.17
0.00
3.010561
-2.420128
1.450000
-4.145597
78
195.08
-2.420128
3.789750
0.770510
0.25
1.15
Al
2.863924
1.403115
20.418205
23.195740
6.613165
3.527021
0.314873
0.365551
0.379846
0.759692
-2.807602
-0.301435
1.258562
-1.247604
-2.83
0.00
0.622245
-2.488244
0.785902
-2.824528
13
26.981539
-2.488244
3.527021
0.759692
0.85
1.15
Pb
3.499723
0.647872
8.450154
8.450063
9.121799
5.212457
0.161219
0.236884
0.250805
0.764955
-1.422370
-0.210107
0.682886
-0.529378
-1.44
0.00
0.702726
-0.538766
0.935380
-1.439436
82
207.2
-0.538766
5.212457
0.764955
0.85
1.15
Fe
2.481987
1.885957
20.041463
20.041463
9.818270
5.236411
0.392811
0.646243
0.170306
0.340613
-2.534992
-0.059605
0.193065
-2.282322
-2.54
0.00
0.200269
-0.148770
0.391750
-2.539945
26
55.847
-0.148770
5.236411
0.340613
0.85
1.15
Mo
2.728100
2.723710
29.354065
29.354065
8.393531
4.476550
0.708787
1.120373
0.137640
0.275280
-3.692913
-0.178812
0.380450
-3.133650
-3.71
0.00
0.875874
0.776222
0.790879
-3.712093
42
95.94
0.776222
4.476550
0.275280
0.85
1.15
Ta
2.860082
3.086341
33.787168
33.787168
8.489528
4.527748
0.611679
1.032101
0.176977
0.353954
-5.103845
-0.405524
1.112997
-3.585325
-5.14
0.00
1.640098
0.221375
0.848843
-5.141526
73
180.9479
0.221375
4.527748
0.353954
0.85
1.15
W
2.740840
3.487340
37.234847
37.234847
8.900114
4.746728
0.882435
1.394592
0.139209
0.278417
-4.946281
-0.148818
0.365057
-4.432406
-4.96
0.00
0.661935
0.348147
0.582714
-4.961306
74
183.84
0.348147
4.746728
0.278417
0.85
1.15
Mg
3.196291
0.544323
7.132600
7.132600
10.228708
5.455311
0.137518
0.225930
0.5
1.0
-0.896473
-0.044291
0.162232
-0.689950
-0.90
0.00
0.122838
-0.226010
0.431425
-0.899702
12
24.305
-0.226010
5.455311
1.0
0.85
1.15
Co
2.505979
1.975299
27.206789
27.206789
8.679625
4.629134
0.421378
0.640107
0.5
1.0
-2.541799
-0.219415
0.733381
-1.589003
-2.56
0.00
0.705845
-0.687140
0.694608
-2.559307
27
58.9332
-0.687140
4.629134
1.0
0.85
1.15
Ti
2.933872
1.863200
25.565138
25.565138
8.775431
4.680230
0.373601
0.570968
0.5
1.0
-3.203773
-0.198262
0.683779
-2.321732
-3.22
0.00
0.608587
-0.750710
0.558572
-3.219176
22
47.88
-0.750710
4.680230
1.0
0.85
1.15
Zr
3.199978
2.230909
30.879991
30.879991
8.559190
4.564902
0.424667
0.640054
0.5
1.0
-4.485793
-0.293129
0.990148
-3.202516
-4.51
0.00
0.928602
-0.981870
0.597133
-4.509025
40
91.224
-0.981870
4.564902
1.0
0.85
1.15
+22 −0
Original line number Diff line number Diff line
EAM database tool
Xiaowang Zhou (Sandia), xzhou at sandia.gov

based on this paper:

X. W. Zhou, R. A. Johnson, and H. N. G. Wadley, Phys. Rev. B, 69,
144113 (2004).

This tool can be used to create an DYNAMO-formatted EAM
setfl file for alloy systems, using any combination 
of the elements discussed in the paper and listed in
the EAM_code file, namely:

Cu, Ag, Au, Ni, Pd, Pt, Al, Pb, Fe, Mo, Ta, W, Mg, Co, Ti, Zr

Steps:

1) compile create.f -> a.out  (e.g. gfortran create.f)
2) edit the input file EAM.input to list 2 or more desired elements to include
3) a.out < EAM.input will create an EAM *.set file
4) in DYNAMO or LAMMPS lingo, this is a setfl file
   that can be used with the LAMMPS pair_style eam/alloy command
+248 −0
Original line number Diff line number Diff line
C author: X. W. Zhou, xzhou@sandia.gov
c      open(unit=5,file='a.i')
      call inter
c      close(5)
      call writeset
      stop
      end
ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
c main subroutine.                                                c
ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
      subroutine inter
      character*80 atomtype,atommatch,outfile,outelem
      namelist /funccard/ atomtype
      common /pass1/ re(16),fe(16),rhoe(16),alpha(16),
     *   beta(16),beta1(16),A(16),B(16),cai(16),ramda(16),
     *   ramda1(16),Fi0(16),Fi1(16),Fi2(16),Fi3(16),
     *   Fm0(16),Fm1(16),Fm2(16),Fm3(16),Fm4(16),
     *   fnn(16),Fn(16),rhoin(16),rhoout(16),rhol(16),
     *   rhoh(16),rhos(16)
      common /pass2/ ielement(16),amass(16),Fr(5000,16),
     *   rhor(5000,16),z2r(5000,16,16),ntypes,blat(16),
     *   nrho,drho,nr,dr,rc,outfile,outelem
      ntypes=0
10    continue
      atomtype='none'
      read(5,funccard)
      if (atomtype .eq. 'none') goto 1200
      open(unit=10,file='EAM_code',form='FORMATTED',status='OLD')
11    read(10,9501,end=1210)atommatch
9501  format(a80)
      if (atomtype .eq. atommatch) then
         ntypes=ntypes+1
         length=len_trim(outfile)
         if (length .eq. len(outfile)) then
            outfile = atomtype
         else
            outfile = outfile(1:length)//atomtype
         endif
         length=len_trim(outelem)
         if (length .eq. len(outelem)) then
            outelem = atomtype
         else
            outelem = outelem(1:length)//' '//atomtype
         endif
         read(10,*) re(ntypes)
         read(10,*) fe(ntypes)
         read(10,*) rhoe(ntypes)
         read(10,*) rhos(ntypes)
         read(10,*) alpha(ntypes)
         read(10,*) beta(ntypes)
         read(10,*) A(ntypes)
         read(10,*) B(ntypes)
         read(10,*) cai(ntypes)
         read(10,*) ramda(ntypes)
         read(10,*) Fi0(ntypes)
         read(10,*) Fi1(ntypes)
         read(10,*) Fi2(ntypes)
         read(10,*) Fi3(ntypes)
         read(10,*) Fm0(ntypes)
         read(10,*) Fm1(ntypes)
         read(10,*) Fm2(ntypes)
         read(10,*) Fm3(ntypes)
         read(10,*) fnn(ntypes)
         read(10,*) Fn(ntypes)
         read(10,*) ielement(ntypes)
         read(10,*) amass(ntypes)
         read(10,*) Fm4(ntypes)
         read(10,*) beta1(ntypes)
         read(10,*) ramda1(ntypes)
         read(10,*) rhol(ntypes)
         read(10,*) rhoh(ntypes)
         blat(ntypes)=sqrt(2.0)*re(ntypes)
         rhoin(ntypes)=rhol(ntypes)*rhoe(ntypes)
         rhoout(ntypes)=rhoh(ntypes)*rhoe(ntypes)
      else
         do 1 i=1,27
1        read(10,*)vtmp
         goto 11
      endif
      close(10)
      goto 10
1210  write(6,*)'error: atom type ',atomtype,' not found'
      stop
1200  continue
      nr=2000
      nrho=2000
      alatmax=blat(1)
      rhoemax=rhoe(1)
      do 2 i=2,ntypes
         if (alatmax .lt. blat(i)) alatmax=blat(i)
         if (rhoemax .lt. rhoe(i)) rhoemax=rhoe(i)
2     continue
      rc=sqrt(10.0)/2.0*alatmax
      rst=0.5
      dr=rc/(nr-1.0)
      fmax=-1.0
      do 3 i1=1,ntypes
      do 3 i2=1,i1
      if ( i1 .eq. i2) then
         do 4 i=1,nr
            r=(i-1.0)*dr
            if (r .lt. rst) r=rst
            call prof(i1,r,fvalue)
            if (fmax .lt. fvalue) fmax=fvalue
            rhor(i,i1)=fvalue
            call pair(i1,i2,r,psi)
            z2r(i,i1,i2)=r*psi
4        continue
      else
         do 5 i=1,nr
            r=(i-1.0)*dr
            if (r .lt. rst) r=rst
            call pair(i1,i2,r,psi)
            z2r(i,i1,i2)=r*psi
            z2r(i,i2,i1)=z2r(i,i1,i2)
5        continue
      endif
3     continue
      rhom=fmax
      if (rhom .lt. 2.0*rhoemax) rhom=2.0*rhoemax
      if (rhom .lt. 100.0) rhom=100.0
      drho=rhom/(nrho-1.0)
      do 6 it=1,ntypes
      do 7 i=1,nrho
         rhoF=(i-1.0)*drho
         call embed(it,rhoF,emb)
         Fr(i,it)=emb
7     continue
6     continue
      return
      end
ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
c This subroutine calculates the electron density.                c
ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
      subroutine prof(it,r,f)
      common /pass1/ re(16),fe(16),rhoe(16),alpha(16),
     *   beta(16),beta1(16),A(16),B(16),cai(16),ramda(16),
     *   ramda1(16),Fi0(16),Fi1(16),Fi2(16),Fi3(16),
     *   Fm0(16),Fm1(16),Fm2(16),Fm3(16),Fm4(16),
     *   fnn(16),Fn(16),rhoin(16),rhoout(16),rhol(16),
     *   rhoh(16),rhos(16)
      f=fe(it)*exp(-beta1(it)*(r/re(it)-1.0))
      f=f/(1.0+(r/re(it)-ramda1(it))**20)
      return
      end
ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
c This subroutine calculates the pair potential.                  c
ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
      subroutine pair(it1,it2,r,psi)
      common /pass1/ re(16),fe(16),rhoe(16),alpha(16),
     *   beta(16),beta1(16),A(16),B(16),cai(16),ramda(16),
     *   ramda1(16),Fi0(16),Fi1(16),Fi2(16),Fi3(16),
     *   Fm0(16),Fm1(16),Fm2(16),Fm3(16),Fm4(16),
     *   fnn(16),Fn(16),rhoin(16),rhoout(16),rhol(16),
     *   rhoh(16),rhos(16)
      if (it1 .eq. it2) then
         psi1=A(it1)*exp(-alpha(it1)*(r/re(it1)-1.0))
         psi1=psi1/(1.0+(r/re(it1)-cai(it1))**20)
         psi2=B(it1)*exp(-beta(it1)*(r/re(it1)-1.0))
         psi2=psi2/(1.0+(r/re(it1)-ramda(it1))**20)
         psi=psi1-psi2
      else
         psi1=A(it1)*exp(-alpha(it1)*(r/re(it1)-1.0))
         psi1=psi1/(1.0+(r/re(it1)-cai(it1))**20)
         psi2=B(it1)*exp(-beta(it1)*(r/re(it1)-1.0))
         psi2=psi2/(1.0+(r/re(it1)-ramda(it1))**20)
         psia=psi1-psi2
         psi1=A(it2)*exp(-alpha(it2)*(r/re(it2)-1.0))
         psi1=psi1/(1.0+(r/re(it2)-cai(it2))**20)
         psi2=B(it2)*exp(-beta(it2)*(r/re(it2)-1.0))
         psi2=psi2/(1.0+(r/re(it2)-ramda(it2))**20)
         psib=psi1-psi2
         call prof(it1,r,f1)
         call prof(it2,r,f2)
         psi=0.5*(f2/f1*psia+f1/f2*psib)
      endif
      return
      end
ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
c This subroutine calculates the embedding energy.                c
ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
      subroutine embed(it,rho,emb)
      common /pass1/ re(16),fe(16),rhoe(16),alpha(16),
     *   beta(16),beta1(16),A(16),B(16),cai(16),ramda(16),
     *   ramda1(16),Fi0(16),Fi1(16),Fi2(16),Fi3(16),
     *   Fm0(16),Fm1(16),Fm2(16),Fm3(16),Fm4(16),
     *   fnn(16),Fn(16),rhoin(16),rhoout(16),rhol(16),
     *   rhoh(16),rhos(16)
      if (rho .lt. rhoe(it)) then
         Fm33=Fm3(it)
      else 
         Fm33=Fm4(it)
      endif
      if (rho .lt. rhoin(it)) then
         emb=Fi0(it)+
     *       Fi1(it)*(rho/rhoin(it)-1.0)+
     *       Fi2(it)*(rho/rhoin(it)-1.0)**2+
     *       Fi3(it)*(rho/rhoin(it)-1.0)**3
      else if (rho .lt. rhoout(it)) then
         emb=Fm0(it)+
     *       Fm1(it)*(rho/rhoe(it)-1.0)+
     *       Fm2(it)*(rho/rhoe(it)-1.0)**2+
     *       Fm33*(rho/rhoe(it)-1.0)**3
      else
         emb=Fn(it)*(1.0-fnn(it)*log(rho/rhos(it)))*
     *       (rho/rhos(it))**fnn(it)
      endif
      return
      end
ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
c write out set file.                                             c
ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
      subroutine writeset
      character*80 outfile,outelem
      common /pass1/ re(16),fe(16),rhoe(16),alpha(16),
     *   beta(16),beta1(16),A(16),B(16),cai(16),ramda(16),
     *   ramda1(16),Fi0(16),Fi1(16),Fi2(16),Fi3(16),
     *   Fm0(16),Fm1(16),Fm2(16),Fm3(16),Fm4(16),
     *   fnn(16),Fn(16),rhoin(16),rhoout(16),rhol(16),
     *   rhoh(16),rhos(16)
      common /pass2/ ielement(16),amass(16),Fr(5000,16),
     *   rhor(5000,16),z2r(5000,16,16),ntypes,blat(16),
     *   nrho,drho,nr,dr,rc,outfile,outelem
      character*80 struc
      struc='fcc'
      outfile = outfile(1:index(outfile,' ')-1)//'.set'
      open(unit=1,file=outfile)
      write(1,*)
      write(1,*)
      write(1,*)
      write(1,8)ntypes,outelem
8     format(i5,' ',a24)
      write(1,9)nrho,drho,nr,dr,rc
9     format(i5,e24.16,i5,2e24.16)
      do 10 i=1,ntypes
      write(1,11)ielement(i),amass(i),blat(i),struc
      write(1,12)(Fr(j,i),j=1,nrho)
      write(1,12)(rhor(j,i),j=1,nr)
10    continue
11    format(i5,2g15.5,a8)
12    format(5e24.16)
      do 13 i1=1,ntypes
      do 13 i2=1,i1
      write(1,12)(z2r(i,i1,i2),i=1,nr)
13    continue
      close(1)
      return
      end
Baidu
map