# Intersecting ellipsoids

First of all I know that finding whether two ellipsoids intersect is complicated and there are specific algorithms created for this purpose (see, for instance, . Actually, I am trying to follow this documentation in order to implement something in Mathematica. Of course there are other approaches).

I have "googled" for a while, searching for Mathematica codes but I could not find anything relevant. Maybe this has to do with the fact that functions like RegionIntersection were introduced rather recently (2014; version 10.0).

Anyway, here is a simple code. It's not panacea of course. (Found here)

RegionIntersectQ[e1_, e2_] :=
RegionQ[DiscretizeRegion[RegionIntersection[e1, e2]] // Quiet]


It works in many cases. E.g.

Graphics3D[{g1 =
Ellipsoid[{0, 0, 0}, {{5, 2, 3}, {2, 3, 2}, {3, 2, 5}}],
g2 = Ellipsoid[{1, 1, 1}, {{10, 2, 3}, {2, 3, 2}, {3, 2, 5}}]},
Axes -> True]
RegionIntersectQ[g1, g2] // Timing
(*2.62...,True*) My first question is about the time efficiency. Why 2.6 sec are needed?

We can find (in version 10.3 at least!) other cases that the code finds correctly whether the two ellipsoids intersect (the code tries harder when there is a tiny region of intersection).

But there are cases where it fails strangely and I do not understand why.

Graphics3D[{{Opacity[0.6],
g1 = Ellipsoid[{0, 0, 0}, {{5, 2, 3}, {2, 3, 2}, {3, 2, 5}}]},
g2 = Ellipsoid[{4.4, 1, 1}, {{10, 2, 3}, {2, 3, 2}, {3, 2, 7}}]},
Axes -> True]
RegionIntersectQ[g1, g2] // Timing I am not expecting other people to spend their time to implement algorithms of others but I am wondering how we can use efficient the current functionality of Mathematica (built-in functions like RegionIntersection) to check if two ellipsoids intersect or not.

EDIT After writing the question I found this which has connection with my previous question here. @JasonB addresses the issue very thoroughly.

But the query remains the same: Can I expect Mathematica's relevant functionality to do the job (even with a bit help of us) or should I forget the whole idea?

• Just checking, but you do know RegionIntersection is not implemented for MeshRegion objects embedded in 3D? – M.R. Dec 4 '15 at 9:05
• @M.R. I don't see your point. Can you elaborate more on your comment? Thanks. – Dimitris Dec 4 '15 at 9:13
• Check the docs for ref/RegionIntersection under Possible issues, and it says RegionIntersection functionality in 3d is not complete. – M.R. Dec 4 '15 at 9:18
• BTW, @dimitris, what are you working on with all these graphics problems? They are really quite fun and I'm tempted to blow off a days work to implement the algorithms you linked to. – Jason B. Dec 4 '15 at 10:13
• @Saurav: I want to determine if there is intersection. No need for analytic expressions. – Dimitris Dec 7 '15 at 14:38

## 2 Answers

FindInstance is faster:

FindInstance[RegionMember[g1, {x, y, z}] && RegionMember[g2, {x, y, z}],
{x, y, z}] // AbsoluteTiming
(*  {0.002422, {{x -> 0, y -> 0, z -> 0}}}  *)

FindInstance[
RegionMember[
Ellipsoid[{5, 0, 0}, {{5, 2, 3}, {2, 3, 2}, {3, 2, 5}}], {x, y, z}] && RegionMember[
Ellipsoid[{1, 1, 1}, {{10, 2, 3}, {2, 3, 2}, {3, 2, 5}}], {x, y, z}],
{x, y, z}] // AbsoluteTiming
(*  {0.045233, {}}  *)


Also, Reduce is faster than the OP's approximate discretizing, but a bit slower than the above in the positive (first) case.

Reduce[Exists[{x, y, z}, RegionMember[g1, {x, y, z}] && RegionMember[g2, {x, y, z}]],
{x, y, z}] // AbsoluteTiming
(*  {0.011632, True}  *)

Reduce[Exists[{x, y, z},
RegionMember[Ellipsoid[{5, 0, 0}, {{5, 2, 3}, {2, 3, 2}, {3, 2, 5}}], {x, y, z}] &&
RegionMember[Ellipsoid[{1, 1, 1}, {{10, 2, 3}, {2, 3, 2}, {3, 2, 5}}], {x, y, z}]],
{x, y, z}] // AbsoluteTiming
(*  {0.044602, False}  *)

• I saw the answer right now. Thank you very much! – Dimitris Jan 14 '17 at 12:15

Starting in M11.1 you can use RegionDisjoint. Here are the two OP examples:

g1 = Ellipsoid[{0, 0, 0}, {{5, 2, 3}, {2, 3, 2}, {3, 2, 5}}];
g2 = Ellipsoid[{1, 1, 1}, {{10, 2, 3}, {2, 3, 2}, {3, 2, 5}}];

RegionDisjoint[g1, g2] //AbsoluteTiming


{0.002053, False}

g1 = Ellipsoid[{0, 0, 0}, {{5, 2, 3}, {2, 3, 2}, {3, 2, 5}}];
g2 = Ellipsoid[{4.4, 1, 1}, {{10, 2, 3}, {2, 3, 2}, {3, 2, 7}}];

RegionDisjoint[g1, g2] //AbsoluteTiming


{0.001416, False}

And here is @Michael's disjoint example:

r1 = Ellipsoid[{5, 0, 0}, {{5, 2, 3}, {2, 3, 2}, {3, 2, 5}}];
r2 = Ellipsoid[{1, 1, 1}, {{10, 2, 3}, {2, 3, 2}, {3, 2, 5}}];

RegionDisjoint[r1, r2] //AbsoluteTiming


{0.041539, True}