Can someone please explain the difference between

fix ID group-ID rigid/small ...
fix ID group-ID nvt temp ...

fix ID group-ID rigid/nvt/small ...
Basically I want to simulate CO2 molecule, but when I use the earlier method it becomes numerically unstable and when I use the later one bond angle changes from 180deg. Is there anything else I should do with fix rigid? (Maybe putting angle_coeff zero? This doesn’t seem right!)
PS: I am using traPPE model of CO2.

Please have a close look at your log file and observe any warnings printed by LAMMPS.
That should give you a crucial hint.

Or to put it more succinctly, one is correct use and one is incorrect.

1 Like

Got it. I used to ignore the incosistent image flags, which was messing with the simulation.

That is not the correct answer to your question.

When using fix nvt and fix rigid on the same atoms, you get a conflict because both fixes apply time integration so that the both try to update atom positions and velocities and that leads to inconsistent data and bogus trajectories and eventually a crash.

LAMMPS prints a warning in these cases.

That either means, your input is not correct or there are other mistakes.

You’re right! I set the image flags for molecule but after a small simualtion, angle’s still changing. Now the problem is my log file doesn’t have any warnings to fix the input script. Maybe the potentials are wrong? But I am using a standard traPPE model. I’m sharing my case file and molecule file along with this, If there is any obvious mistake, can you please check? Otherwise, I’ll debug it one line at a time.
Thanks a lot for your help!
Testrun.rar (19.8 KB)

I think I found the problem. I am using angle_style harmonic with angle_coeff 0 and 180. I think the first parameter in this coeff is harmonic ‘spring’ constant?

If you are using fix rigid correctly, bond and angle potentials should not matter.

I just looked at your input and you are not using fix rigid correctly. Your input does not keep any molecules rigid until after your fix npt run. Up until that point your CO2 molecules can bend and stretch which they are not supposed to do for the TraPPE model, if I remember correctly.
This also causes problems with minimization, because rigid fixes are not invoked during minimization, so the best you can do it short runs with regular quenching of the kinetic energy and damped dynamics (i.e. using fix viscous to slow high energy atoms down).

1 Like

Thank you for your suggestions. I fixed my case using fix deposit to create molecules and then fix viscous reduced the KE of the system. I am getting good trajectories now!
I just have one doubt. When we perform equilibration, which parameters should we look for? Since my system is quite large, I’m getting quite high total energy, but the fluctuation are almost zero in it and I’m using fix NVT for equilibration, so the density remains constant. That’s why I’m not sure, if my system is actually equilibrated. (Thermostat is also quite stable +/-10K)

Technically speaking, anything and everything. In practice, you can never know whether your system is trapped in some (deep) well in the potential hyper-surface and whether there may be an even deeper well elsewhere, only that it is difficult to reach. Therefore it is very system specific and it requires significant experience to make a choice when to say that the risk of missing out is low.
This is where you need to consult with somebody with significant experience in your specific area of research. It is not something where one can just give a “do this, not that” kind of rigid and general advice.

1 Like