#### Minimal Enclosing Circle

```/*
* Minimal Enclosing Circle
*   The minimal enclosing circle of a set of points is the circle with the
*   smallest radius containing all the points
*/

/* mec(Points, p(Xc,Yc), RadiusSquared) is true if the circle whose centre   */
/*   is at (Xc,Yc) and whose radius squared is RadiusSquared is the minimum  */
/*   enclosing circle of the set of Points.                                  */
/*   This is Emo Welzl's O(n) time algorithm.                                */
/* e.g. mec([p(3,1),p(4,8),p(3,10),p(8,8)], C, R) gives C=p(4.1,5.5),R=21.46 */
mec([Point], Point, 0.0):-!.

/* mec_1(OtherPoints, BoundaryPoints, Centre, RadiusSquared) is true if the  */
/*   BoundaryPoints are on the circle with the given Centre and              */
/*   RadiusSquared, and the OtherPoints are on or inside this same circle.   */
circle_3(R1, R2, R3, Centre, RadiusSquared), !.
mec_1(_, [_,_,_], p(0,0), 0.0):-!,
write(`Two identical points`), nl,
abort.
distance_squared(P, Centre, DistanceSquared),

/* circle_2(p(X1,Y1), p(X2,Y2), p(Xc,Yc), RadiusSquared) is true if the      */
/*   centre of the smallest circle passing through the two points (X1,Y1)    */
/*   and (X2,Y2) is (Xc,Yc), and RadiusSquared is the square of its radius.  */
/* e.g. circle_2(p(-6,8), p(-2,-4), p(X,Y), R) gives X=-4, Y=2, R=40         */
Xc is (X1+X2) / 2.0,
Yc is (Y1+Y2) / 2.0,
X1_X2 is X1 - X2,
Y1_Y2 is Y1 - Y2,
RadiusSquared is (X1_X2*X1_X2 + Y1_Y2*Y1_Y2) / 4.0.

/* circle_3(p(X1,Y1), p(X2,Y2), p(X3,Y3), p(Xc,Yc), RadiusSquared) is true   */
/*   if the centre of circle passing through the three points (X1,Y1),       */
/*   (X2,Y2) and (X3,Y3) is (Xc,Yc), and RadiusSquared is the square of its  */
/*   radius. Fails if the three points are collinear.                        */
/*   This algorithm is borrowed from the comp.graphics.algorithms FAQ.       */
/* e.g. circle_3(p(-6,5), p(-3,-4), p(2,1), p(X,Y), R) gives X=-3, Y=1, R=25 */
A is X2 - X1,
B is Y2 - Y1,
G is 2.0 * (A*(Y3 - Y2) - B*(X3 - X2)),
G =\= 0.0,
C is X3 - X1,
D is Y3 - Y1,
E is A*(X1 + X2) + B*(Y1 + Y2),
F is C*(X1 + X3) + D*(Y1 + Y3),
Xc is (D*E - B*F) / G,
Yc is (A*F - C*E) / G,
H is X1 - Xc,
I is Y1 - Yc,