Support Sympy for tabulated interactions

[RudolfWeeber]

It would be very handy, if one cut set up tabulated interactions sympy expressions. Using Sympy, it is also possible to obtain the force automatically.

import sympy as sp

r = sp.Symbol('r)

h_energy = \
  epsilon *sp.Piecewise(((1-r/sigma)**(5/2),r<sig"),(0,True))
bond = espressomd.interactions.tabulated_from_sympy(
    type="distance", var=r, energy=h_energy,min=0,max=10,steps=1000)

Implementation

• use sp.diff() w.r.t var to obtain force from energy
• use sp.lambdify() to generate a callable function
• use np.linspace to generate the distances between min and max
• calculate force and energy arrays
• put them into the desired tabulated bond

Example implementation

d = sp.Symbol('d')

h_energy = \
  force_model["hertzian_eps"] *sp.Piecewise(((1-d/force_model["hertzian_sig"])**(5/2),d<force_model["hertzian_sig"]),(0,True))
h_force = -sp.diff(h_energy,d)
h_energy = sp.lambdify(d, h_energy,"numpy")
h_force = sp.lambdify(d,h_force,"numpy")
print(f"{h_energy=}, {h_force=}")
dists = np.linspace(0,system.box_l[0],5000)
hertz_bond = espressomd.interactions.TabulatedDistance(
   min=0, max=system.box_l[0], energy=h_energy(dists), force=h_force(dists))