Dear all,
I have been trying to implement Ewald method to alpha-quartz. I referred to this two posts:
http://lammps.sandia.gov/threads/msg13329.html
The main idea is to get relaxed structure of crystal alpha-quartz at 0K, quenching from ambient environment. Basically, I used tabulated BKS + repulsive potential for short range (<10.17) with PPPM=1.0e-6 and set g-ewald parameter in k_space modify. In my understanding, g-ewald = alpha. Now I want to set R_cut=3.2/alpha=10.17 while I do not know how to set R_cut in script. Does it mean R_cut is automatically chosen to be the outer-cutoff (5.5) in tabulated potential?
About reciprocal vector k=2pi(n_x/L_x, n_y/L_y, n_z/L_z), I want to sum up to n_max=alphaL_x, alphaL_y, alpha*L_z separately in three direction since my simulation region is a cuboid. I have not found a way to set kmax yet and even log.lammps did not tell what kmax they took.
Thanks
Bingyu Cui
Input scripts:
#Initialization############################################################
units metal
dimension 3
boundary p p p
atom_style charge
Atom definition##########################################################
lattice custom 5.4047 &
a1 0.9092 0.0000 0.0000 &
a2 -0.4546 0.7873 0.0000 &
a3 0.0000 0.0000 1.0000 &
basis 0.4697 0.0000 0.0000 &
basis 0.0000 0.4697 0.6667 &
basis 0.5303 0.5303 0.3333 &
basis 0.4133 0.2672 0.1188 &
basis 0.2672 0.4133 0.5479 &
basis 0.7328 0.1461 0.7855 &
basis 0.5867 0.8539 0.2145 &
basis 0.8539 0.5867 0.4521 &
basis 0.1461 0.7328 0.8812
region simbox block 0.0 10.0 0.0 10.0 0.0 10.0 units lattice
create_box 2 simbox
create_atoms 1 box basis 1 1 basis 2 1 basis 3 1 basis 4 2 basis 5 2 basis 6 2 basis 7 2 basis 8 2 basis 9 2
mass 1 28.0855
mass 2 15.9994
group siliconatoms type 1
group oxygenatoms type 2
set group siliconatoms charge 2.4
set group oxygenatoms charge -1.2
Atoms interactions settings##################################
Si type 1, O type 2
pair_style table linear 30262 pppm
pair_coeff 1 1 potential.Ewald SiSi
pair_coeff 1 2 potential.Ewald SiO
pair_coeff 2 2 potential.Ewald OO
kspace_style pppm 1.0e-6
kspace_modify gewald 0.31465093412
neighbor 2.0 bin
neigh_modify every 1 delay 0 check yes
pair_write 1 2 10001 r 0.001 5.5 tableEW.txt Si_O 2.4 -1.2
pair_write 2 2 10001 r 0.001 5.5 tableEW.txt O_O -1.2 -1.2
pair_write 1 1 10001 r 0.001 5.5 tableEW.txt Si_Si 2.4 2.4
Equilibrate
velocity all create 298.0 635845 dist gaussian
fix 0 all npt temp 298.0 0.1 0.1 iso 1.01325 0.0 1.0
thermo_style custom step temp pe etotal press vol density lx ly lz pxx pyy pzz pxy pxz pyz fmax
thermo 10000
dump 1 all atom 100000 dump1.qua
run 100000
fix 4 all box/relax x 0.0 y 0.0 z 0.0 vmax 0.005
min_style cg
minimize 1e-25 1e-25 5000 10000
Potential:
SiSi
N 30262 RSQ 0.001 5.501
1 0.001 -2.310806983947287e19 -1.3864848852556973e23
2 0.031638584039112745 -2.2502271011385094e10 -4.3606647477092646e12
3 0.044732538492690085 -2.4459180019968233e9 -3.801479245700687e11
…
30260 5.500818211866304 0. 0.
30261 5.5009091066840945 0. 0.
30262 5.501 0. 0.
SiO
N 30262 RSQ 0.001 5.501
1 0.001 1.4187525355067962e36 1.7025030426081557e40
2 0.031638584039112745 1.4102697348252677e18 5.348923579191188e20
3 0.044732538492690085 2.2101620569234056e16 5.929005054657236e18
…
30261 5.5009091066840945 0. 0
30262 5.501 0. 0
OO
N 30262 RSQ 0.001 5.501
1 0.001 1.4200972364567373e38 1.7041166837480846e42
2 0.031638584039112745 1.4116063957330952e20 5.3539933164696e22
…
30261 5.5009091066840945 0. 0
30262 5.501 0. 0