question about pair_coul_long_soft potential

You are mentioning the Coulomb part of pair soft. What is calculated in the pair styles is
(-r du/dr) / r^2
so the exponent of the denominator is in this case 3/2. Since denc = sqrt(…) was needed before, that’s where the cube comes from.

Maybe doing things like this is old fashioned and should be left for compilers to optimise. If you can do faster, please contribute!