Smallest enclosing circle problem: Difference between revisions

</pre>

=={{header|Phix}}==

Based on the same code as Wren, and likewise limited to 2D circles - I believe (but cannot prove) the main barrier to coping with more than two dimensions lies wholly within the circle_from3() routine.

<lang Phix>type point(sequence) return true end type

enum CENTRE, RADIUS

type circle(sequence) return true end type

function distance(point a, b)

return sqrt(sum(sq_power(sq_sub(a,b),2)))

end function

function in_circle(point p, circle c)

return distance(p,c[CENTRE]) <= c[RADIUS]

end function

function circle_from2(point a, b)

-- return the smallest circle that intersects 2 points:

point centre = sq_div(sq_add(a,b),2) -- (ie midpoint)

atom radius = distance(a,b)/2

circle res = { centre, radius }

return res

end function

function circle_from3(point a, b, c)

-- return a unique circle that intersects three points

atom {{aX,aY},{bX,bY},{cX,cY}} = {a,sq_sub(b,a),sq_sub(c,a)},

B = bX * bX + bY * bY,

C = cX * cX + cY * cY,

D = bX * cY - bY * cX

point centre = { (cY * B - bY * C) / (2 * D) + aX,

(bX * C - cX * B) / (2 * D) + aY }

circle res = { centre, distance(centre, a) }

return res

end function

function trivial(sequence r)

integer l = length(r)

switch l do

case 0: return {{0,0},0}

case 1: return {r[1],0}

case 2: return circle_from2(r[1],r[2])

case 3: return circle_from3(r[1],r[2],r[3])

end switch

end function

function welzlr(sequence p, r)

if p={} or length(r)=3 then return trivial(r) end if

point p1 = p[1]

p = p[2..$]

circle d = welzlr(p, r)

if in_circle(p1,d) then return d end if

return welzlr(p, append(r,p1))

end function

procedure welzl(sequence p)

printf(1,"Points %v ==> centre %v radius %.14g\n",prepend(welzlr(shuffle(p),{}),p))

end procedure

constant tests = {{ { 0, 0 }, { 0, 1 }, { 1, 0 } },

{ { 5, -2 }, { -3, -2 }, { -2, 5 }, { 1, 6 }, { 0, 2 } }}

papply(tests,welzl)</lang>

{{out}}

<pre>

Points {{0,0},{0,1},{1,0}} ==> centre {0.5,0.5} radius 0.70710678118655

Points {{5,-2},{-3,-2},{-2,5},{1,6},{0,2}} ==> centre {1,1} radius 5