# Ellipsoid¶

`Elllipsoids`

are defined by their local coordinates and three `semi-axes`

.

## Ellipsoid-Ellipsoid contact¶

Thsi algorithm is implemented in by `woo.dem.Cg2_Ellipsoid_Ellipsoid_L6Geom`

.

Following [PRPraestgaardL96], we note `semi-axes`

\(a_i\) and \(b_i\) as \(i\in\{1,2,3\}\); local axes are given as sets of orthonormal vectors \(\vec{u}_i\), \(\vec{v}_i\). The intercenter vector reads \(\vec{R}=\vec{r}_b-\vec{r}_a\), with \(\vec{r}_a\) and \(\vec{r}_b\) being positions of ellipsoid centroids. We define matrices

which are both invertible for non-vanishing semi-axes; then we define the function

where \(\lambda\in[0,1]\) is an unknown parameter. The function \(S(\lambda)\) is concave and non-zero for \(\lambda\in[0,1]\) and zero at both end-points. It therefore has a sigle maximum. Using the Brent’s method, we obtain the maximum value \(\Lambda\) for which \(S(\Lambda)\) is maximal.

Todo

This computation can be (perhaps substantially) sped up by using other iterative methods, such as Newton-Raphson, finding root of \(S'(\lambda)\), and re-using the value from the previous step as the initial guess. There are also papers suggesting better algorithms such as [ZIPM09].

We define the Perram-Wertham potential (first introduced in [PW85]) as

for which ([PW85])

\(F_{AB}<1\) if both ellipsoids overlap,

\(F_{AB}=1\) if they are externaly tangent,

\(F_{AB}>1\) if the ellipsoids do not overlap.

[PW85] gives a geometrical interpretation of \(F_{AB}\), where

where \(\mu\) is scaling factor which must be applied to both ellipses to be tangent.

Following [DTS05] eq (18), we can compute the contact normal and the contact point of two ellipsoids as

with \(\mat{G}\) defined in (1); note that \(\vec{n}_c\) is not normalized.

The penetration depth (overlap distance) can be reasoned out as follows. \(\mu\) scaled ellipsoid sizes while keeping their distance, so that they become externally tangent. Therefore \(1/\mu\) scales ellipsoid distance while keeping their sizes. With \(d=|\mat{R}|\) being the current inter-center distance, we obtain

This is the displacement that msut be performed along \(\vec{R}\) while the contact normal may be oriented differently; we therefore project \(u_n'\) along \(\vec{R}\) onto (normalized) \(\vec{n}_c\) obtaining

The \(u_n\), \(\normalized{\vec{n}}\), \(\vec{r}_c\) can be fed to Generic contact routine for further computation.