special_bonds command?

Dear all

I am simulating a system including TIP3P water model. In this system I use fix shake to constrain hydrogen atoms as well as H-O-H angles ( … m 1.008 a 1). I also use :neigh_modify delay 2 every 1.

I would appreciate if you let me know if I have to use special_bonds command?
If yes, what should be its keywords and values.

Thanks a lot for your attention,

Hi Sally
   Fix shake should work independently of the special_bonds settings.
Any choice of special_bonds settings (including the default) should
not effect fix shake's behavior.
   Of course, if you are NOT using fix shake, then the special_bonds
setting would be important. (In that case, the two H atoms might
exert a repulsive force on each other. If you water model was not
parameterized to take this into account, then you would need to use
something like "special bonds lj/coul 0 0 0" or "special bonds lj/coul
0 0 1" to insure that you turn off the pair interactions between the H
atoms.) But if you ARE successfully using fix shake, then you don't
have to worry about this.

    Although I admit that I have sometimes had trouble getting fix
shake to work as I intended. If it helps, here are the values I am
using in the SPCE water model:

    bond_style harmonic
    angle_style harmonic
    bond_coeff 1 1000.0 1.0 # k_bond=1000 does not matter
    angle_coeff 1 1000.0 109.47 # k_angle=1000 does not matter
    group spce type 1 2
    fix fShakeSPCE spce shake 0.0001 10 100 b 1 a 1



Hi Andrew

Thanks a lot for your answer.

I use fix shake successfully for TIP3P water model, In case of F3C
model you know it is flexible and we should not use fix shake, but
according to the following note in LAMMPS

IMPORTANT NOTE: When both a bond and pair style is defined, the
special_bonds command often needs to be used to turn off (or weight)
the pairwise interaction that would otherwise exist between 2 bonded

So, I think even in case of F3C model the default values ( 0 0 0)
should be fine. Am I right?

Thanks a lot,

I think so. If your water model does not use fix shake, then you do
need to think about the special_bonds settings. However I'm not
familiar with the F3C model, but after taking a quick look at it, it
certainly seems as though you should turn off the "1-3" interactions,
either by using the default special_bonds settings (equivalent to
"special_bonds lj/coul 0.0 0.0 0.0"), or "special_bonds lj/coul 0 0
1", as you suggested).



(Incidentally, if you have other molecules in your system which force
you to use non-standard special_bonds settings, OR if there is some
other reason which prevents you from using the "special_bonds"
command, then you can get around this limitation by defining an extra
bond between the two hydrogen atoms of strength 0 (if 1-2 interactions
are disabled), or using a tabular bond between the two H atoms to
cancel the effect of the pair force if you can't disable it. I doubt
this is an issue for most people. Hope this helps.)