FindMinimum and the minimum-radius circle
- To: mathgroup at smc.vnet.net
- Subject: [mg50167] FindMinimum and the minimum-radius circle
- From: DrBob <drbob at bigfoot.com>
- Date: Wed, 18 Aug 2004 01:20:21 -0400 (EDT)
- Reply-to: drbob at bigfoot.com
- Sender: owner-wri-mathgroup at wolfram.com
As you say, a global method shouldn't be necessary; but FindMinimum evidently can't do the job.
Reducing the number of starting values for each variable to one AND specifying Method->QuasiNewton fails miserably in my tests, even though I'm using the convex hull and even with only three data points. I always get a FindMinimum::lstol error message (totally spurious, I think), and I often get a circle that can't be optimum.
I tried specifying all the available methods, and they're all equally bad.
Needs["Statistics`"]
Needs["Graphics`"]
Needs["DiscreteMath`ComputationalGeometry`"]
sq = #1 . #1 & ;
sqDiff[x_, y_] = sq[{x, y} - #1] & ;
sqDiff[{x_, y_}] = sq[{x, y} - #1] & ;
diff = Abs[(#1 - #2)/(#1 + #2)] < 0.0001 & ;
circleFinder[n_Integer] :=
Module[{data, hull, r, pt, x2, y2, radius},
data = RandomArray[NormalDistribution[0, 1],
{n, 2}]; hull = data[[ConvexHull[data]]];
radius[x_, y_] = Max[sqDiff[x, y] /@ hull];
{x2, y2} = Median[hull]; {r, pt} =
FindMinimum[radius[x, y], {{x, x2}, {y, y2}},
Method -> ConjugateGradient]; r = Sqrt[r];
pt = {x, y} /. pt; {data, hull, r, pt,
Length[hull] + 1 - Length[
Union[sqDiff[pt] /@ hull, SameTest -> diff]]}]
plotter[n_] := Module[{data, hull, r, pt, count},
{data, hull, r, pt, count} = circleFinder[n];
Print[count]; Show[Graphics[{PointSize[0.02],
Point /@ data, Red, Point[pt], Circle[pt, r],
Blue, Line[Join[hull, {First[hull]}]],
Point /@ hull}], AspectRatio -> Automatic]]
counter[n_] := Last[circleFinder[n]]
circleFinder[3]
NMinimize with constraints is the clear winner so far.
Bobby
Thomas E Burton <tburton at brahea.com> wrote in message news:<cfsi5h$9st$1 at smc.vnet.net>...
> (See thread "Smalest enclosing circle".)
>When I read Bobby's updated message changing FindMinimum to NMinimize,I thought, "There must be only one local minimum. Why do we need aglobal method?". Indeed, a contour plot of the function radius[x,y]shows exactly one minimum. Because the smallest circle in most casestouches three points, the contours surrounding the minimum aretriangular. FindMinimum as implemented by Bobby (with two startingvalues for each coordinate) gets stuck at the first acute corner ithits. (And maybe some obtuse corners as well--I have not thought thisthrough.)
>Though frequently maligned in this group for being too brief, theImplementation Notes do shed light on this issue: When two startingvalues are given for each variable, FindMinimum defaults to Brent'sprincipal axis method. If you specify Method->"QuasiNewton", or simplyprovide only one starting value for each coordinate instead of two,then FindMinimum does a respectable job in a small fraction of the timeneeded for NMinimize. On the other hand, if you follow Bobby's trackand reduce the field of points to its convex hull, chances are that youan easily afford either method.
>Tom Burton