2
$\begingroup$

I have two dynamic perpendicular parabolas: $y=ax^2+c, x=by^2+d$, which intersect in four points (the coefficients are chosen that way). There is a fact that we can construct a circle through these four points.

My question is which command allows to construct a circle through four points?

$\endgroup$

closed as off-topic by m_goldberg, MarcoB, user9660, Bob Hanlon, Yves Klett May 17 '16 at 4:35

This question appears to be off-topic. The users who voted to close gave this specific reason:

  • "This question cannot be answered without additional information. Questions on problems in code must describe the specific problem and include valid code to reproduce it. Any data used for programming examples should be embedded in the question or code to generate the (fake) data must be included." – m_goldberg, MarcoB, Community, Bob Hanlon, Yves Klett
If this question can be reworded to fit the rules in the help center, please edit the question.

  • 1
    $\begingroup$ Three points determine a circle. $\endgroup$ – Mark Adler May 15 '16 at 20:00
  • $\begingroup$ yes, and i can graph it. But because of the machine accuracy being nonezero it can pass near the fourth point $\endgroup$ – Michael Freimann May 15 '16 at 20:19
  • 1
    $\begingroup$ @Michael Freimann If you edit your post to include the smallest and nicest (indented four spaces to let SE know this is code and not to be desktop published) bit of code, NOT a screenshot image, that shows your two parabolas with sliders then someone else might take a few minutes to try to write a couple of lines of code that would try to fit the best circle onto your four points. $\endgroup$ – Bill May 15 '16 at 20:48
  • 1
    $\begingroup$ "There is a fact that we can construct a circle through these four points." - there's a constructive proof of that statement; why not use Mathematica to work through it? $\endgroup$ – J. M. is away May 16 '16 at 2:29
17
$\begingroup$

This is a motivating post. I have made no effort to deal with intersection <4 pts. I post only illustrative examples.

The centre of the circle is the intersection of perpendicular bisectors of chords.

f[a_, b_, c_, d_] := 
 Module[{e1 = a x^2 + b, e2 = c y^2 + d, s, p, l, ln, ctr},
  s = Solve[{x, e1} == {e2, y}, {x, y}, Reals];
  p = {#, a #^2 + b} & /@ N[x /. s];
  l = Partition[p, 2, 1];
  ln = InfiniteLine[{(#1 + #2)/2, (#1 + #2)/2 + Cross[#2 - #1]}] & @@@
     l;
  ctr = RegionIntersection @@ ln[[{1, 2}]];
  ParametricPlot[{{x, e1}, {e2, y}}, {x, -3, 3}, {y, -3, 3}, 
   Epilog -> {Red, PointSize[0.03], Point[p], Black, Dashed, ln, 
     Green, ctr, Purple, Thick, 
     Circle[ctr[[1]], Norm[ctr[[1]] - p[[1]]]]}, 
   AspectRatio -> Automatic, 
   PlotLabel -> Row[{"{", e1, ",", e2, "}"}]]]

Examples:

Grid[Partition[
  f @@@ {{1, -2, 1, -2}, {1, -3, -1, 3}, {1, -4, 1, -2}, {-1, 2, 
     1, -3}}, 2]]

enter image description here

I reiterate I have made no effort to deal with number of intersections <4.

$\endgroup$
  • $\begingroup$ ...and of course, Circumsphere[] is built-in, and works in 2D. :) $\endgroup$ – J. M. is away May 16 '16 at 9:31
  • $\begingroup$ @J.M. Thanks...didn't know that...always learning. :) $\endgroup$ – ubpdqn May 16 '16 at 9:34
  • $\begingroup$ It's a convenient function, but I think your approach is more pedagogically appropriate for the OP. $\endgroup$ – J. M. is away May 16 '16 at 9:45
6
$\begingroup$

Here is a standard algebraic approach, to complement @ubpdqn's geometric one.

If we have a sequence of points $(x_i,y_i)$, $i=1,2,\dots,n$ lying on a circle, then the linear system $$A\,(x_i^2+y_i^2)+B\,x_i+C\,y_i+D=0,\quad i=1,2,\dots,n$$ will have a nontrivial solution $(A_0,B_0,C_0,D_0)$. If there are at least three distinct points not on a line, the solution space will have dimension 1 and the circle through the points will be unique.

Here is a function that will find the intersections of curves and see a circle runs through the points.

ClearAll[circleIntersect];
circleIntersect::nocirc = "There is not a unique circle through the points.";
circleIntersect[eqs_, vars : {_, _}] := 
 Module[{circle, mat, form, res},
  form = {vars.vars, Sequence @@ vars, 1};       (* {x^2+y^2, x, y, 1} *)
  mat = form /. Solve[eqs, vars, Reals];         (* NSolve is often much faster *)
  Check[                   (* Dot::dotsh, Set::shape errors if null space not 1D *)
   {circle} = Sow[NullSpace[mat], "basis"].form; (* make basis available via Reap *)
   res = circle == 0,
   Message[circleIntersect::nocirc];
   res = $Failed
   ];
  Sow[mat[[All, {2, 3}]], "points"];             (* make points available via Reap *)
  res /; FreeQ[res, $Failed]
  ]

Examples

Two parabolas intersecting in four points.

p1 = y == x^2 - 4;
p2 = x == 2 y^2 - 3;
{eq, {pts}} = Reap[circleIntersect[{p1, p2}, {x, y}], "points"];
If[FreeQ[eq, circleIntersect],
 ContourPlot[Evaluate@{p1, p2, eq}, {x, -4, 4}, {y, -4, 4},
  Epilog -> {Red, PointSize[Medium], Point@pts}]]

Mathematica graphics

Solve will return points in terms of Root objects. If there is a nice form for the solution, FullSimplify might find it:

eq // FullSimplify
(*  2 (x^2 + y^2) == 11 + x + 2 y  *)

Tangent parabolas require some care. In this case if the machine real 4. is replaced by the exact integer 4, Solve runs for longer than I'm willing to wait. Numeric error with the machine-precision 4. leads to nonzero imaginary parts in the solution; hence the use of Chop below. One could replace 4. by the arbitrary-precision 4.`16, which would lead to a zero imaginary part.

p1 = y == x^2 - 4.;
p2 = x == y^2 - Root[16411 - 1152 #1 - 4096 #1^2 + 256 #1^3 &, 2];
{eq, {pts}} = Reap[circleIntersect[{p1, p2}, {x, y}], "points"];
If[FreeQ[eq, circleIntersect],
 ContourPlot[Evaluate@{p1, p2, Chop@eq}, {x, -4, 4}, {y, -4, 4},
  Epilog -> {Red, PointSize[Medium], Point@pts}]]

Mathematica graphics

eq
(*
  (0.960573 + 9.29229*10^-16 I) + (0.16052 + 1.32742*10^-14 I) x +
     (0.16052 - 1.06369*10^-14 I) y - (0.16052 - 6.13072*10^-15 I) (x^2 + y^2) == 0
*)

Points not lying on a circle. This intersection of a parabola and a quartic does not define a circle.

p1 = y == x^2 - 4 - x^4/16;
p2 = x == 2 y^2 - 3;
circleIntersect[{p1, p2}, {x, y}]
If[FreeQ[%, circleIntersect],
 ContourPlot[Evaluate@{p1, p2, %}, {x, -4, 4}, {y, -4, 4}]]

Mathematica graphics

More than one circle. This intersection of two parabolas has only two points and does not define a unique circle.

p1 = y == x^2 - 4;
p2 = x == 2 y^2 - 1;
{eq, {{{pts}}, {{basis}}}} = 
  Reap[circleIntersect[{p1, p2}, {x, y}], {"points", "basis"}];
pts // N

Mathematica graphics

Note, one could detect this. The basis returned NullSpace consists of two vectors. Linear combinations of these parameterize all the circles through the two points of intersection.

Show[
 ContourPlot[Evaluate@{p1, p2}, {x, -4, 4}, {y, -4, 4}, ContourStyle -> Black],
 ContourPlot[  (* family of circles *)
  Evaluate[
   Table[{(1 - t^2)/(1 + t^2), (2 t)/(1 + t^2)}, {t, 1.8^Range[-8, 8]}] . 
    N@basis.(List /@ {x^2 + y^2, x, y, 1})], 
  {x, -4, 4}, {y, -4, 4},
  Contours -> {{0}}]
 ]

Mathematica graphics

$\endgroup$
  • $\begingroup$ beautiful...+1 of course :) $\endgroup$ – ubpdqn May 17 '16 at 2:47

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