Hi

Dear Agraj Abhishek,

Please include this mail in lammps forum so that everyone can be benefit.

See below for the input script which you requested. I take the simulation box boundary lx (from thermo_style) as the lattice constant. I wrote a simple code to read lammps output and calculate the lattice constants by averaging over large equilibrium time steps. I believe you could do it too. I don’t understand what you mean by U shape curve, maybe you can post your query to lammps forum.

Good luck,
Christopher

 # strontium titanate system

 atom_style    charge
 units         real
dimension     3
 boundary    p p p

# lattice
lattice        custom 3.9051 a1 1.0 0.0 0.0 a2 0.0 1.0 0.0 a3 0.0 0.0 1.0 basis 0.0 0.0 0.0 basis 0.5 0.5 0.5 basis 0.5 0.5 0.0 basis 0.5 0.0 0.5 basis 0.0 0.5 0.5
region        simbox block 0 9 0 9 0 9
create_box    3 simbox
create_atoms    1 box basis 1 1 basis 2 2 basis 3 3 basis 4 3 basis 5 3
group        sr type 1
group        ti type 2
 group        o type 3
 mass        1 87.62
 mass        2 47.867
mass        3 15.999
set        group sr charge 1.2
set        group ti charge 2.4
set         group o charge -1.2

 # force-field
 pair_style    hybrid/overlay born 11.0 coul/long 11.0 morse 11.0
 pair_coeff     * * coul/long
 pair_coeff     1 1 born 0.32 0.32 2.396 100.0 0.0
 pair_coeff    1 2 born 0.34 0.34 2.253 250.0 0.0
pair_coeff    1 3 born 0.32 0.32 3.124 200.0 0.0
pair_coeff    2 2 born 0.36 0.36 2.110 625.0 0.0
 pair_coeff    2 3 born 0.34 0.34 2.981 500.0 0.0
 pair_coeff    3 3 born 0.32 0.32 3.852 400.0 0.0
 pair_coeff    1 3 morse 2.41 1.18 2.7615
pair_coeff     2 3 morse 4.23 3.82 2.1923
 kspace_style    ewald 0.0001
 pair_modify    tail yes shift no

 # initialize
 timestep    1.0
 velocity    all create 298.0 1234567 mom yes rot yes dist gaussian
 neighbor     4.0 bin
 neigh_modify    every 1 delay 0 check yes

# output
 variable        N equal step
 variable        pote equal pe
 variable        Etotal equal etotal
 variable        T equal temp
 variable        Press equal press
 variable        V equal vol
 variable        kine equal ke
 variable        Lx equal lx/9
 variable        Ly equal ly/9
 variable        Lz equal lz/9
fix             extra all print 1 "${N} ${T} ${V} ${Press} ${kine} ${pote} ${Etotal}" file strontium_titanate.out
fix             cell all print 1 "${N} ${T} ${Lx} ${Ly} ${Lz}" file strontium_titanate.cell
dump            snapshot all custom 1000 strontium_titanate.lammpstrj id q type x y z

 # equilibration & thermalization

fix         NPT all npt temp 298.0 298.0 100.0 iso 1.01325 1.01325 1000.0
run 220000

Thanks for the mail Christopher, I will run the simulation to see the result. Now I have one more doubt, if you can please clarify it, you have done thermalization and equilabration but not time averaging. Is it not required in this case. ? Is it that after so long time all the parameter saturates to a constant value.

Agraj

Hi Agraj,

I don’t use time/ave command to do time averaging. Rather, I write a code to calculate a quantity by averaging over equilibrium time steps, this is called time averaging. From my script, 20 000 is the time steps for the system to reach equilibrium and following 200 000 is the time steps where my system equilibrates.

Regards,
Christopher