Hello Dr. Kohlmeyer,
I have been working with the Li-Al-O system. My simulation consists of minimization followed by NPT runs. The temperature in NPT (T_start = T_stop = 300) is stable at around 300 towards the end of the simulation but the pressure (Pstart = Pstop = 0) is not reaching zero. It is either a large positive or large negative (-7000 to +6000) value even during a long period of time. Am I missing something in the code?
units metal
dimension 3
boundary p p p
atom_style charge
read_data final2.data
variable mLi equal 6.941
variable mAl equal 26.9815
variable mO equal 15.9994
mass 1 {mLi}
mass 2 {mAl}
mass 3 ${mO}
group Li type 1
group Al type 2
group O type 3
set type 1 charge 0.7
set type 2 charge 1.5
set type 3 charge -1.1
pair_style hybrid/overlay table linear 5000 coul/long 10.0
pair_coeff 1 1 table LiAlO2_TOT2009_WithoutCoulomb.table LiLi
pair_coeff 1 2 table LiAlO2_TOT2009_WithoutCoulomb.table LiAl
pair_coeff 1 3 table LiAlO2_TOT2009_WithoutCoulomb.table LiO
pair_coeff 2 2 table LiAlO2_TOT2009_WithoutCoulomb.table AlAl
pair_coeff 2 3 table LiAlO2_TOT2009_WithoutCoulomb.table AlO
pair_coeff 3 3 table LiAlO2_TOT2009_WithoutCoulomb.table OO
pair_coeff * * coul/long
kspace_style pppm 1.0e-4
timestep 0.0001
compute pe all pe/atom
neighbor 1.0 bin
neigh_modify every 5 delay 0 check no
thermo 20
thermo_style custom step dt time temp pe lx ly lz pxx pyy pzz
thermo_modify norm no
min_style cg
minimize 0 0 1000 10000
write_dump all custom dump.end.min id type x y z modify sort id
fix 1 all npt temp 300 300 1 iso 0 0 1 drag 1
thermo 100
thermo_style custom step temp press ke pe etotal vol density
run 1000000
write_dump all custom dump.end id type x y z modify sort id
Sincerely,
Ankit
1 Like
This topic has been discussed on this mailing list many, many times. Please study the mailing list archives.
That is the wrong method. There can be a large hysteresis, even for macroscopic experiments. There are discussions on this topic and using the coexistence method also in the mailing list archives.
Thanks Dr. Kohlmeyer.
I read a few of the mailing lists regarding this problem and learnt that one can check if the average pressure in the last few runs remains zero. I was able to get a zero average pressure taking the average of the last 1000 runs with 0.001 ps timestep.
But currently I am facing another problem while trying to calculate melting temperature via cooling and heating. When I ramp down the temperature from 2500 to 300 K, my system temperature is overshooting to an absurd value of 20000-25000. Is there anything I can rectify in this section?
units metal
dimension 3
boundary p p p
atom_style charge
#atom_modify map hash
read_data final2.data
variable myT equal 573
variable mLi equal 6.941
variable mAl equal 26.9815
variable mO equal 15.9994
mass 1 {mLi}
mass 2 {mAl}
mass 3 ${mO}
group Li type 1
group Al type 2
group O type 3
set type 1 charge 1.75
set type 2 charge 5.25
set type 3 charge -3.5
#Li-Li
pair_style hybrid/overlay buck/coul/long 13 5 zbl 0.5 2
pair_coeff 1 1 buck/coul/long 0 0.285710 0
pair_coeff 1 1 zbl 3 3
Yes I read about hysteresis in this method. But I just wanted to get an approximate range of the melting point.
The problem I am facing is that the temperature is overshooting to 25K or 30K while ramping down from 2500 to 300 Kelvin. Is it because of the potential? The density that I obtained using this potential was very close to the experimental value.
If your model is even somewhat representative of real physics (and that’s not guaranteed), then the process of freezing your system should release a large amount of latent heat (enthalpy of fusion).
How much heat is released, and over what time period? Should you expect your thermostat to be able to remove the heat that quickly? What happens to your system if it can’t?
Do you know of a prior paper which uses a similar method to determine the melting point of a similar system?
Cheers,
Shern
Add-on question: is the force field parameterized to correctly model the liquid phase?
Dr. Kohlmeyer and Dr. Tee,
So it turned out that the force field was not appropriate for the system as the radial distribution function of the solid phase from this force field was not representative of a crystal. I have adopted a different force field taken from “Atomistic simulations of the defect chemistry and self-diffusion of Li-ion in LiAlO2.” Energies 12, no. 15 (2019): 2895. " and now I am using the solid liquid interface method to find out the melting point. Thank you for your valuable inputs.