10
$\begingroup$

I have a parametric equation of two ellipses as follow

$\begin{cases} x=a_1 \text{sin$\theta $}+b_1 \text{cos$\theta $}+c_1 \\ y=d_1 \text{sin$\theta $}+e_1 \text{cos$\theta $}+f_1 \\ \end{cases}$

$\begin{cases} x=a_2 \text{sin$\theta $}+b_2\text{cos$\theta $}+c_2 \\ y=d_2 \text{sin$\theta $}+e_2 \text{cos$\theta $}+f_2 \\ \end{cases}$

Now I want to solve the intersection point coordinates and visualize them.

My trial

Here, mat1 owns the style {{a1,b1,c1},{d1,e1,f1}}.

solvePoints0[mat1_, mat2_] :=
 Module[{sol, θValue, θ, ptRule},
  sol =
   Solve[
    Thread[
     mat1.{Sin[θ1], Cos[θ1], 1} ==
      mat2.{Sin[θ2], Cos[θ2], 1}], {θ1, θ2}];
  θValue = Mod[(#[[1, 2, 1]] & /@ sol) /. C[1] -> 1., 2 Pi];
  ptRule = List /@ Thread[θ -> θValue];
  Cases[
   mat1.{Sin[θ], Cos[θ], 1} /. ptRule,{_Real, _Real}]
 ]

showEllipse0[mat1_, mat2_, opts : OptionsPattern[]] :=
 With[{pts = solvePoints0[mat1, mat2]},
   ParametricPlot[
    {mat1.{Sin[θ], Cos[θ], 1}, 
     mat2.{Sin[θ], Cos[θ], 1}}, {θ, 0, 2 Pi},
    Epilog -> {PointSize[Large], Red, Point[pts]},
    Evaluate[
     Sequence @@ FilterRules[{opts}, Options[ParametricPlot]]]
   ]
  ] /; MatrixQ[mat1, NumericQ] && MatrixQ[mat2, NumericQ]

Test

showEllipse0[##, ImageSize -> 200] & @@@
 {{{{1, 2, 1}, {2, 3, 4}}, {{2, 2, 2}, {3, 2, 4}}},
  {{{3, 1, 1}, {2, 3, 4}}, {{4, 2, 2}, {3, 4, 4}}},
  {{{3, 2, 1}, {2, 3, 4}}, {{2, 2, 2}, {3, 2, 4}}},
  {{{2, 1, 1}, {2, 3, 4}}, {{4, 1, 2}, {3, 2, 4}}}}

enter image description here

However, solvePoints0 has a bug

solvePoints0[{{4, 2, 1}, {2, 5, 6}}, {{6, 2, 2}, {3, 2, 4}}]
solvePoints0[{{5, 2, 1}, {2, 5, 6}}, {{6, 2, 2}, {3, 2, 4}}]
{}

{}

In fact, they have intersection point coordinates.

showEllipse0 @@@
 {{{{4, 2, 1}, {2, 5, 6}}, {{6, 2, 2}, {3, 2, 4}}},
  {{{5, 2, 1}, {2, 5, 6}}, {{6, 2, 2}, {3, 2, 4}}}}

enter image description here

To fixed this bug, I rewrite the function solvePoints0 using NSolve

solvePoints[mat1_, mat2_] :=
 Module[
   {sol, θValue, θ, ptRule},
   sol =
    Quiet@
      NSolve[
        Thread[
         mat1.{Sin[θ1], Cos[θ1], 1.} ==
          mat2.{Sin[θ2], Cos[θ2], 1.}], {θ1, θ2}] /. C[1] -> 1.;
   θValue = #[[1, 2]] & /@ sol;
   ptRule = List /@ Thread[θ -> θValue];
   Cases[
    mat1.{Sin[θ], Cos[θ], 1} /. ptRule, {_Real, _Real}]
 ] /; MatrixQ[mat1, NumericQ] && MatrixQ[mat2, NumericQ]

Now,

solvePoints[{{4, 2, 1}, {2, 5, 6}}, {{6, 2, 2}, {3, 2, 4}}]
solvePoints[{{5, 2, 1}, {2, 5, 6}}, {{6, 2, 2}, {3, 2, 4}}]

enter image description here

Question

  • I think my methods is fussy and bad, so I would like to know how to solve this question in a elegant way.

Update

Thanks for J.M's suggestions, with the help of Graphics`Mesh`FindIntersections

newSolution[mat1_, mat2_] :=
 Module[{graph, pts},
  graph =
   ParametricPlot[
    {mat1.{Sin[θ], Cos[θ], 1},
     mat2.{Sin[θ], Cos[θ], 1}}, {θ, 0, 2 Pi},Epilog -> {Point[{.1, .2}]}];
  pts = Graphics`Mesh`FindIntersections[First[graph]];
  graph /.
   (Epilog -> _) -> Epilog -> {Red, PointSize[Medium], Point[pts]}
]

Test

newSolution[##] & @@@ 
 {{{{1, 2, 1}, {2, 3, 4}}, {{2, 2, 2}, {3, 2, 4}}}, 
  {{{3, 1, 1}, {2, 3, 4}}, {{4, 2, 2}, {3, 4, 4}}}, 
  {{{3, 2, 1}, {2, 3, 4}}, {{2, 2, 2}, {3, 2, 4}}}, 
  {{{2, 1, 1}, {2, 3, 4}}, {{4, 1, 2}, {3, 2, 4}}}}

enter image description here

newSolution @@@ 
 {{{{4, 2, 1}, {2, 5, 6}}, {{6, 2, 2}, {3, 2, 4}}},
  {{{5, 2, 1}, {2, 5, 6}}, {{6, 2, 2}, {3, 2, 4}}}}

enter image description here

Obviously, this is not a good method owning to that it cannot find the intersection points correctly.

$\endgroup$
  • 6
    $\begingroup$ A sketch: use GroebnerBasis[] to produce the implicit Cartesian equations of the two ellipses, and feed those equations to Solve[]. Retain only the real solutions, and you're done. $\endgroup$ – J. M. is away Jun 25 '15 at 8:03
  • $\begingroup$ Alternatively, if you do not need high accuracy, see this answer. $\endgroup$ – J. M. is away Jun 25 '15 at 8:06
  • 1
    $\begingroup$ Your code cannot figure out what is real due to phantom imaginary parts arising in numerical approximations involving arctans and the like. $\endgroup$ – Daniel Lichtblau Jun 25 '15 at 15:27
  • $\begingroup$ @J. M. Thanks for your link. However, the FindIntersections cannot find all the intersection points correctly(Namely,it has wrong points or lacks of some right points). Please see my update. $\endgroup$ – xyz Jun 26 '15 at 1:56
  • 3
    $\begingroup$ Shutao, note that the wrong points are precisely the points where the ellipse closes up, so you know where they are ($\theta=0$), and you can remove them from the results (with some tolerance since the results are not exactly the same). $\endgroup$ – J. M. is away Jun 26 '15 at 2:16
14
$\begingroup$

Alright, I managed to borrow a computer. Here's an implementation of my suggestion:

ellipseIntersections[mat1_?MatrixQ, mat2_?MatrixQ] /; 
       Dimensions[mat1] == Dimensions[mat2] == {2, 3} :=
{\[FormalX], \[FormalY]} /. 
RootReduce[Solve[Flatten[Map[
      GroebnerBasis[Append[Thread[{\[FormalX], \[FormalY]} == #],
                           \[FormalC]^2 + \[FormalS]^2 == 1],
                    {\[FormalX], \[FormalY]}, {\[FormalC], \[FormalS]}] &,
    {mat1, mat2}.{\[FormalC], \[FormalS], 1}]] == 0,
{\[FormalX], \[FormalY]}, Reals]]

where, apart from the use of GroebnerBasis[], I use formal symbols for safety, and add validity checks for the arguments. I use the combination RootReduce[Solve[(* stuff *)]] to generate exact coordinates; replace this with NSolve[] if you only want numerical approximations.

Before demonstrating this function, allow me to show a different way to render your ellipses. Here, I have chosen to use the NURBS (that is, BSplineCurve[]) representation of a circle, along with a suitable affine transformation. Thus:

makeEllipse[mat_?MatrixQ] /; Dimensions[mat] == {2, 3} := Module[{aff},
    aff = AffineTransform[{Drop[mat, None, -1], mat[[All, -1]]}];
    BSplineCurve[aff /@ {{1, 0}, {1, 1}, {-1, 1}, {-1, 0},
                         {-1, -1}, {1, -1}, {1, 0}},
                 SplineClosed -> True, SplineDegree -> 2, 
                 SplineKnots -> {0, 0, 0, 1/4, 1/2, 1/2, 3/4, 1, 1, 1}, 
                 SplineWeights -> {1, 1/2, 1/2, 1, 1/2, 1/2, 1}]]

Now, pictures!

e1 = {{5, 2, 1}, {2, 5, 6}}; e2 = {{6, 2, 2}, {3, 2, 4}};
Graphics[{{RGBColor[7/19, 37/73, 22/31], makeEllipse[e1]},
          {RGBColor[59/67, 11/18, 1/7], makeEllipse[e2]},
          {Red, AbsolutePointSize[5], Point[ellipseIntersections[e1, e2]]}},
         Frame -> True]

intersecting ellipses

Here's the case e1 = {{1, 2, 1}, {2, 3, 4}}; e2 = {{2, 2, 2}, {3, 2, 4}};: tangent ellipses

which shows that ellipseIntersections can deal with tangencies.

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.