Hello,
I am trying to model a new pair potential in LAMMPS, but I am kind of facing a tricky situation. As far as I have seen, we require the pairwise force and energy expressions for calculations. Here I want to implement a potential which has an energy expression of the form:
E_i = a*sum(f(r_ij)) + b/sum(g(r_ij))
where E_i is the energy associated with the i’th atom and the sum() is over index of neighbours j. It becomes a little tricky to find the derivative (pairwise force) due to the presence of a summation term in the denominator. I did find the derivative:
F(R_ij) = a*f’(R_ij) - (b/sum(g(r_ij))^2)*g’(R_ij)
But implementing the above did not seem to give the force values I am expecting. I guess there might be some issues with double counting, which would have been easy to take care of if there had not been a denominator summation term.
I would like to have some thoughts and suggestions on how I can go about implementing the above function in LAMMPS.
Thanks,
Adarsh