[lammps-users] fix rdf question

I am attempting to output the radial distribution function for a simulation and am using the fix rdf command. I would like to only output the intermolecular rdf and exclude the counting of any atoms within the same molecule. The documentation states that I can avoid counting any bonded atoms (1-2, 1-3, 1-4 interactions within a molecular topology) by using the command “special_bonds 0.0 0.0 0.0”. However, when I do this, the rdf’s are still counting the bonded atoms and are washing out all of the intermolecular data that I am interested in.

Is anyone familiar with using this fix command? Is my interpretation of the “special_bonds” command correct in that if you set the values to 0.0, LAMMPS will not include those bonded atoms when creating the rdf? Are there any other ways, such as post-processing, that have been implemented for computing rdf’s?

I am using the 10Nov05 LAMMPS version for this work. I have attached my input file for reference.

Thanks,

Tim

Initialization

units real

dimension 3

boundary p p p

atom_style full

Hi Tim. Are you using long-range electrostatics? If
so, even with "special_bonds 0.0 0.0 0.0", I think
that the 1-2, 1-3, 1-4 interactions would still be in
the neighbor list (and hence the rdf) and the
documentation should be mention this.

Even so, LAMMPS's fix_rdf is intended for computing
atomic rather than molecular rdfs, so I'd recommend
dumping atom positions and using a post-processor to
compute you're molecule-based rdfs. One way to do this
would be to write a module for pizza.py.

Paul

Paul is correct that using long-range Coulombics
will cause 1-2,1-3,1-4 pairs to be in the neighbor
list that would otherwise not be, and thus mess up
the logic in fix rdf. But it would be an easy addition
to test for this special case and exclude their contribution
to the RDF if the weight = 0.0. Is this the right thing
to do? Does it make sense to count a full RDF contribution
if the weight is less than 1.0 but non-zero?

Steve

My vote would be to turn it off for weight = 0.0 and
on otherwise. But it is somewhat arbitrary --- just
depends what the user wants in the rdf. If the users
want something really complex, it would probably be
better in a post-processor. This fix_rdf was really
only intended to be a very simple rdf calculator for
on-the-fly use. -Paul