How to find the incircle and circumcircle for an irregular polygon?

With the help of J.M. and Chip Hurst, I was able to get this up and running with a reasonably fast (meaning it's still slow) and clean code. The problem of the circumcircle is solved by the new-in-version-10.4 function BoundingRegion, and the incircle is found by using NMinimize along with SignedRegionDistance.

In the code below, I overload the definitions to account for the fact that when you call CountryData[<country>,"Polygon"] you get your coordinates wrapped in GeoPosition, where the x and y coordinates are backwards. But if you call CountryData[<country>,{"Polygon","Mercator"}] then it's a normal polygon.

I also wrote the code so that it will find the total circumcircle for a collection of polygons, so it takes account for all the land masses in a disjointed geographic region.

ClearAll[inCircle, outCircle];
inCircle[Polygon[GeoPosition[a_]]] := 
  inCircle[Polygon[a /. {x_?NumericQ, y_?NumericQ} :> {y, x}]];
inCircle[Polygon[a_ /; Depth[a] == 4]] := 
  Last@(SortBy[inCircle /@ a, Last]);
inCircle[a_List] := inCircle[Polygon@a];
inCircle[Polygon[a_ /; Depth[a] == 3]] := inCircle[Polygon[a]] =
  Module[{interior},
    interior = BoundaryDiscretizeGraphics@Polygon[a];
    {{x, y} /. #2, Abs@#1} & @@ 
     NMinimize[
      SignedRegionDistance[interior, {x, y}], {x, y} ∈ 
       interior]
    ];

outCircle[p_List] := List @@ BoundingRegion[p, "MinDisk"];
outCircle[Polygon[GeoPosition[a_]]] := 
  outCircle[Polygon[a /. {x_?NumericQ, y_?NumericQ} :> {y, x}]];
outCircle[Polygon[a_]] := outCircle[Flatten[a, Depth[a] - 3]];

Just for educational value, I also made a version of this that runs on Mathematica 9, posted as a separate answer. Here is the result for Japan,

Show[Graphics[{Blue, Disk @@ outCircle[#]}], Graphics@#, 
   Graphics[{Red, Disk @@ inCircle[#]}]] &@
 CountryData["Japan", {"Polygon", "Mercator"}]

And here I sort the countries by their roundness:

countries = CountryData[#, "StandardName"] & /@ CountryData[];
resultsList = Table[
      With[{pol = CountryData[country, {"Polygon", "Mercator"}]},
       {CountryData[country, "Name"], Last[inCircle[pol]]^2/
        Last[outCircle[pol]]^2, 
        Show[Graphics[{Blue, Disk @@ outCircle[pol]}], Graphics@pol, 
         Graphics[{Red, Disk @@ inCircle[pol]}]]}],
      {country, countries}] // SortBy[#[[2]] &] // 
    Reverse;~Monitor~country
resultsList[[;;20]] // TableForm

The 20 roundest countries are below, and the full list is posted on gist, and can be downloaded via << "https://git.io/v6mt2"

For the in-circle aspect of the question, we can use the fact that the center of the in-circle is the Voronoi node with the most negative signed distance to the polygon. Finding the Voronoi mesh and the negative signed distance of all of the Voronoi nodes are very fast:

incircle[poly_List?MatrixQ] := With[
    {nodes = VoronoiMesh[poly]["Coordinates"], sd = SignedRegionDistance[Polygon[poly]]},
    With[
        {ord = First @ Ordering[sd[nodes], 1]},

        Circle[
            nodes[[ord]],
            Abs @ sd[nodes[[ord]]]
        ]
    ]
]

Let's compare timings with the accepted answer:

poly = Entity["Country", "France"]["Polygon"][[1,1,1]];

inCircle[poly] //AbsoluteTiming
incircle[poly] //AbsoluteTiming

{2.41236, {{46.7076, 2.52144}, 3.47992}}

{0.116144, Circle[{46.7076, 2.52146}, 3.47992]}

Let me present another approach. The code should be self-explanatory.

country = "Mauritius";
imgwidth = 500; (* increase this for better fit of inner circle *)
allpols = Polygon /@ First @ CountryData[country, {"Polygon", "Mercator"}];
areas = Area /@ allpols;
pol = First @ Pick[allpols, areas, Max[areas]]; (* just the biggest connected land *)
polpoints = First @ pol;
img = Composition[
  ImageCrop,
  Image[#, ImageSize -> imgwidth] &,
  Graphics
] @ pol
skeleton = Composition[
  ImageData,
  ImageAdjust,
  SkeletonTransform,
  ColorNegate
] @ img;
Image @ skeleton (* just to show the skeleton *)
centre = FirstPosition[skeleton, 1.] (* in 'List coordinates' *)
dim = Dimensions[skeleton]
newcentre = {Last @ centre, First @ dim - First @ centre} (* in Cartesian coordinates *)
{xmin, xmax} = Through[{Min, Max}[First /@ polpoints]]
ymin = Min[Last /@ polpoints]
scalefactor = First @ ImageDimensions[img] / First @ Differences[{xmin, xmax}]
scaler = ScalingTransform[ConstantArray[scalefactor, 2], {0, 0}]
translationvector = - {xmin, ymin}
translator = TranslationTransform[translationvector]
newpolpoints = (scaler @* translator) /@ polpoints;
boundary = MeshRegion[
  newpolpoints,
  Line /@ Partition[Append[Range[Length @ newpolpoints], 1], 2, 1]
];
innerradius = RegionDistance[boundary, newcentre]
outercircle = BoundingRegion[newpolpoints, "MinDisk"]
outerradius = Last @ outercircle
Graphics[{
  Blue,
  outercircle,
  White,
  Polygon[newpolpoints],
  Red,
  Disk[newcentre, innerradius]
}]
roundness = (innerradius/outerradius)^2;
Print["Roundness of ", country, ": ", roundness]

Notes

1. pol is taken to be the biggest connected land of a country, in the selected projection. For Japan in the Mercator projection, for example, pol is Hokkaido (which is actually smaller than Honshu). You may want to choose a representative area differently.

2. For some countries like Fiji, you would really need to crank up imgwidth.