I am trying to understand and implement the Gay-Berne potential, but I do not quite understand how one goes about constructing them. I am following this paper to parse through the idea.
I understand how Lennard-Jones/Coulombic potentials work in MD engines: every particle has some $\epsilon, \sigma, q$ associated with it, and these potentials are purely a function of the distance $r=|\mathbf{r}_{ij}|$ between two particles. It does not depend upon the orientation of the two particles. $$U = U(r)$$
In the above paper, a particle $i$ at position $\mathbf{r}_i$ and orientation $\omega$ with respect to the laboratory frame is represented by a Gaussian function $G_i (\omega _i, \mathbf{r})$ which is $$ G_i(\omega _i, \mathbf{r}) = \exp [-\frac{1}{2}(\mathbf{r}-\mathbf{r}_i)^T M_i ^T S_i ^{-2} M_i (\mathbf{r}-\mathbf{r}_i)]$$ where S is a diagonal 'shape' matrix with elements $\sigma _x, \sigma _y, \sigma _z$, the axes of the ellipsoid representing the molecule measures with respect to a unit length $\sigma _0$, and $M_i$ is the rotation matrix transforming from laboratory to a molecular frame.
What I don't understand is, why do we need a Gaussian function to denote a particle? A 3-dimensional ellipse can be represented with the shape matrix S, and its orientation which can be given by the matrix M. Where does $\omega$ figure into the above definition? If anybody has any recommendations on how to go about understanding the implementation of this idea, I would greatly appreciate it.