Hi all.
Is there any reason why it would be difficult to report eangle on a per-atom basis in dump custom? ebond is reported; why not eangle? The energy would be associated with the central atom about which the angle is defined.
--Craig
Craig,
There is no reason other than that your assumption that only the central
atom contributes to the atoms eangle is incorrect. An angle energy
calculation is in essence somewhat of a three-body potential, albeit a very
localized one. This means that 1/3 of the energy calculated for one item in
the angle list should go to each of the contributing atoms. To calculate
the full contribution to one atom, you will have to find the two other
contributing items in the angle list that have the atom at hand in it. You
would have to create an eangle array on the Atom object, which you can
populate by adding 1/3 of the energy for each item calculated.
Pieter
Since angle::single() routines were added in a recent patch, this
would be possible. I.e. write a compute_eangle class like
compute_ebond. This is a short term solution since eventually
we'll accumulate per-atom energy for all interactions including
many-body ones.
I agree with Pieter that the energy of
one angle should be split across 3 atoms, like it is across 2 for
a bond.
Steve
Hi Pieter.
Craig,
There is no reason other than that your assumption that only the central
atom contributes to the atoms eangle is incorrect. An angle energy
calculation is in essence somewhat of a three-body potential, albeit a very
localized one. This means that 1/3 of the energy calculated for one item in
the angle list should go to each of the contributing atoms.
Seems like that is a convention which makes equally as much sense as putting all of the energy on the central atom. (It's not like pair-wise interactions where there is symmetry between the 2 players... in a triplet-forming an angle, the central guy is clearly special). I just thought it would be easier housekeeping to do the latter. If you want to split it 3 ways, that's fine by me. I was just wondering why it's not already built-in. Seems useful.
To calculate
the full contribution to one atom, you will have to find the two other
contributing items in the angle list that have the atom at hand in it. You
would have to create an eangle array on the Atom object, which you can
populate by adding 1/3 of the energy for each item calculated.
Presumably I could just mimic what's done for ebond in ComputeEbondAtom::compute_peratom(). I'd need to write a "compute_eangle_atom.cpp". It seems like I would allocate a data array and associate it with Compute::scalar_atom like in compute_ebond_atom.cpp. Does that sound right?
--Craig
Craig,
I'm afraid it is a little more than a convention. It has to do with how
forces and energies should be thought of when distributing them on a per
atom basis. You wouldn't put all of the energy calculated for a pair on
only one of the contributing atoms either.
Your approach on how to add it to a dump file sounds like a reasonable
approach.
Pieter