Greetings,
I’d like to write a new angle style that has some similarities to the class2 angle style. I’ve been looking at the code from “angles_class2.cpp” in order to get an idea of how I should write the new angle style. I have a few questions about how the forces are tabulated for the bond angle term E_ba. Specifically, I don’t understand why the term ‘aa11’ is defined as:
aa11 = aa1 * c/rsq1
and not
aa11 = aa1*(1/r1)
Likewise, I don’t understand why ‘aa12’ is defined as:
aa12 = -aa1 / (r1 * r2)
Actually, I don’t even understand what the ‘aa12’ is doing there at all- I would have thought that the ‘aa11’ term would be sufficient to differentiate the theta term in E_ba w.r.t. x . It looks as if there’s an issue related to dealing with the specific symmetry of this particular potential that I’ve neglected. Anyone out there who has experience with the code for this angle style want to answer this? Apologies if this question has already been answered on the mail list, I searched but didn’t find anything.
Thanks
PS This code comes from version LAMMPS-28May12