Dear LAMMPS Users,
I am trying to replicate the system described in http://dx.doi.org/10.1103/PhysRevB.86.045213 however I am having trouble obtaining the correct elastic constants.
It is PbTe with NaCl-type structure.
Potential includes Coulomb part, Morse part between Pb-Te and Te-Te pairs, and cosine/squared interaction for Te-Pb-Te angle.
pair_style hybrid/overlay morse/smooth/linear 8 coul/cut 20.0
pair_coeff * * coul/cut 20.0
pair_coeff 1 2 morse/smooth/linear 0.465 0.863 3.68 7
pair_coeff 2 2 morse/smooth/linear 0.494 1.51 4.22 8
angle_coeff 1 0.340 90.0
I am setting the system via read_restart command from custom made file (PbTe with NaCl-type structure and angles: twelve 90 degree angles Te-Pb-Te per every Pb atom).
When I try to calculate the elastic constants, the minimization (fix 3 all box/relax aniso 0.0) does not reach pressure close to 0 (i.e. it is usually lower than -100 000 bars and does not depend on initial lattice parameter), despite zeroing forces. Also, the minimization itself sometimes converges to unphysical energy values (e.g. 1e+200 eV).
I am using the ELASTIC example.
I have tried different cutoff distances, initial lattice parameters, types of minimization (cg, sd) however without success.
I would be very gratefull if someone could tell me the possible reasons for such behaviour.