1

I have two circles in 3D space generated from their parametric equations. Also, these two circles are on the same plane. How can I find out their intersection points from the parametric equations?

Thanks everyone, your answers are really helpful. However, I just realize I may post this on the wrong place. I am using Matlab not Mathematica. Meanwhile as required, I edit my question to make it more clear. First, these two circles are on the same plane and they have two intersection points. second, the radius, the centers and the normal vector are known. The parametric equation I use to generate the circle is:

$$\quad \quad P(t) = r\cos(t)\,u+r\sin(t) n \otimes u+C$$

where, $C$ is the center, $r$ is the radius, $n$ is the normal vector, and $\otimes$ represents the cross product.

For my purposes, $n = (0.7071, 0.7071, 0.3400), $ $u = (n_y, -n_z, 0)$. The centers are $C1 = (-382.6075, 531.1854, 203.9778), C2 = (-382.6075, 522.1854, 203.9778)$.

The radius is $50$ for both circles.

I will try to work out the matlab code based on your reply.

Asaf Karagila
  • 393,674
96742zs
  • 65

3 Answers3

10

Since you haven't provided any details, I'll provide some examples:

Here are parametric equations of two circles in the same plane in 3D that intersect:

cir1 = {20 Cos[t] + 15, 20 Sin[t], 2};
cir2 = {30 Cos[t], 30 Sin[t], 2};

Here they are visually:

pl = ParametricPlot3D[{cir1, cir2}, {t, -Pi, Pi}, BoxRatios -> {1, 1, 1}]

Mathematica graphics

Let's find their points of intersection. We create parametric regions:

reg1 = ParametricRegion[cir1, {{t, -Pi, Pi}}];
reg2 = ParametricRegion[cir2, {{t, -Pi, Pi}}];

Solve for the points:

sol = Reduce[{x, y, z} ∈ reg1 && {x, y, z} ∈ reg2, {x, y, z}]

Mathematica graphics

We can extract these points:

pts = {x, y, z} /. {ToRules @ sol}

OR

pts = {x, y, z} /. Solve[{x, y, z} ∈ reg1 && {x, y, z} ∈ reg2, {x, y, z}]

Mathematica graphics

Visualize:

Show[pl, Graphics3D[{PointSize[0.02], Red, Point[pts]}]]

Mathematica graphics

  • +1 of course - But for the Solve-variant you should use f.e. ... // Values. Otherwise pts won't plot. –  Sep 18 '14 at 11:48
  • @eldo. Thanks, I forgot about that. I'll update. –  Sep 18 '14 at 11:49
5

Just because this is a different approach, I'm adding another answer. As before here are the circles:

cir1 = {20 Cos[t] + 15, 20 Sin[t], 2};
cir2 = {30 Cos[t], 30 Sin[t], 2};

And the parametric regions:

reg1 = ParametricRegion[cir1, {{t, -Pi, Pi}}];
reg2 = ParametricRegion[cir2, {{t, -Pi, Pi}}];

Here are the points of intersection:

pts = MeshCoordinates @ DiscretizeRegion @ RegionIntersection[reg1, reg2]
{{24.1666667, -17.7756075, 2.}, {24.1666667, 17.7756075, 2.}}

Visualize:

Show[ParametricPlot3D[{cir1, cir2}, {t, -Pi, Pi}, BoxRatios -> {1, 1, 1}], 
 Graphics3D[{PointSize[0.02], Red, Point[pts]}]]

Mathematica graphics

1

worth noting this problem is readily worked out in closed form. Here I'm assuming the circles are coplanar, you would obviously need to verify that if its not known a priori.

randomly rotate the plane of the circles. Dealing with the resulting floating point parametrization makes the problem slightly more challenging...

 trans = RotationTransform[RandomReal[{0, 2 Pi}] , 
           RandomReal[{-1, 1}, 3], {0, 0, 0}];
 cir1 = trans@{20 Cos[t] + 15, 20 Sin[t], 2};
 cir2 = trans@{30 Cos[t], 30 Sin[t], 2};

 cen1 = ((cir1 /. t -> 0) + (cir1 /. t -> Pi))/2;
 cen2 = ((cir2 /. t -> 0) + (cir2 /. t -> Pi))/2;
 cc = Norm[cen1 - cen2];
 r1 = Norm[ (cir1 /. t -> 0) - (cir1 /. t -> Pi) ]/2;
 r2 = Norm[ (cir2 /. t -> 0) - (cir2 /. t -> Pi) ]/2;
 (* check here if not Abs[r1-r2]<=cc<=r1+r2 no intersect *)
 u1 = (cen2 - cen1)/cc;
 n1 = Cross[(cir1 /. t -> 0) - cen1, (cir1 /. t -> Pi/2) - cen1]/r1^2
 u2 = Cross[u1, n1];
 x1 = ( r1^2 - r2^2 + cc^2) /2 /cc;
 x2 = Sqrt[ r1^2 - x1^2 ];
 Show[ParametricPlot3D[{cir1, cir2}, {t, -Pi, Pi}, 
        BoxRatios -> {1, 1, 1}], Graphics3D[{Red, PointSize[.02],
          Point[cen1 + x1 u1 + x2 u2  ] ,
          Point[cen1 + x1 u1 - x2 u2  ] }]]

enter image description here

Edit: graphic illustration of the approach: note by construction u1,u2,n1 is an orthogonal coordinate system because the circles are known to be coplanar. Additional checks need to be done if the circles are not known to be coplanar ( u1.n1 == 0 for example )

enter image description here

george
  • 131
  • 5
  • thanks for the reply, is it possible to explain the math behind it or provide a link ? I am not a mathematica user at all and trying to make a matlab code based on your post. From what I understand, the cen1 and cen2 are the centers of the circle after RotationTransform. cc is the normal vector between two centers. I don't know what the '(cir1 /. t -> 0)' mean. In my case, these two circles have the same radius and their centers and normal vector are known. – 96742zs Sep 19 '14 at 06:41
  • added an illustration that should help explain things. See the docs for basic syntax issues, http://reference.wolfram.com/language/ref/ReplaceAll.html – george Sep 19 '14 at 14:07