Largest Nonoverlapping Smoke Cloud
National Olympiad — National round (2020 season)
Submission name: incinerator.cpp
Problem summary
Existing nonoverlapping circles have centers inside the 1000×1000 square. Place a new circle with its center inside that square, allowing its cloud outside the square, and maximize its radius without interior overlap.
Input
N, followed by center x, center y and radius for each existing circle.
Output
The maximum radius, with absolute error below 10⁻⁹.
Constraints
2 ≤ N ≤ 100; 0 ≤ center coordinates ≤ 1000; 1 ≤ radius ≤ 1414; input decimals have at most 9 fractional digits.
Worked example
An independently prepared example for this explanation.
Input
7 42.1337 13.666 191.17 543.21 765.432 55.555 987.987 666.666 321.321 0 1000 100 666.666 66.66 16.6 133.331 666.666666666 123.321 888.8 88.8 8.8
Output
320.337196456704
Solution analysis
- For a proposed new radius r, inflate each old circle by r. A feasible center is a point of the square outside every inflated disk. Thus feasibility is monotone in r and can be binary-searched.
- If feasible points exist, one exists at a square corner, an intersection of a circle boundary with a square side, or an intersection of two circle boundaries. To see this, follow the boundary of a feasible component: a boundary arc either meets another boundary or is a complete circle. In the complete-circle case it encloses the forbidden disk rather than an isolated feasible component, and any feasible region outside it either reaches the square boundary or another circle.
- Enumerate all these candidate points for each radius and test each against all inflated circles. Use long double, a small geometric tolerance, and ninety binary-search iterations. In the example the best center is near (434.939950686,341.298204535), with radius about 320.337196456704.
Why it works
Nonoverlap is equivalent to distance from every old center being at least oldRadius+r, exactly the inflated-disk complement. The arrangement of circle boundaries and square edges subdivides the square into cells of constant coverage. A nonempty uncovered cell has a closure vertex unless its boundary consists only of complete disjoint circles; in that exceptional configuration its uncovered exterior reaches a square corner or a boundary intersection. Therefore the candidate test detects feasibility. Monotone binary search converges to the largest feasible radius; long-double calculations and the stated tolerance keep numerical slack well below 10⁻⁹ at this coordinate scale.
Watch out. The new cloud may leave the country; constrain its center only. Use long double throughout. Circle tangencies and intersections on the square boundary need inclusive comparisons.
Commented C++17 solution
An independently written complete program with standard headers.
#include<algorithm>
#include<cmath>
#include<iomanip>
#include<iostream>
#include<vector>
// Algorithm: For a proposed new radius r, inflate each old circle by r.
using namespace std;
struct Circle{
long double x,y,radius;
};
int main(){
int n;
cin>>n;
vector<Circle>circles(n);
for(Circle&circle:circles){
cin>>circle.x>>circle.y>>circle.radius;
}
auto feasible=[&](long double radius){
auto clear=[&](long double x,long double y){
const long double tolerance=1e-12L;
if(x<-tolerance||x>1000+tolerance||y<-tolerance||y>1000+tolerance){
return false;
}
for(Circle circle:circles){
long double required=circle.radius+radius-tolerance;
long double dx=x-circle.x;
long double dy=y-circle.y;
if(dx*dx+dy*dy<required*required){
return false;
}
}
return true;
};
for(int x:{0,1000}){
for(int y:{0,1000}){
if(clear(x,y)){
return true;
}
}
}
for(int i=0;i<n;i++){
Circle a=circles[i];
long double ra=a.radius+radius;
for(int side:{0,1000}){
long double square=ra*ra-(side-a.x)*(side-a.x);
if(square>=0){
long double offset=sqrtl(square);
if(clear(side,a.y+offset)||clear(side,a.y-offset)){
return true;
}
}
square=ra*ra-(side-a.y)*(side-a.y);
if(square>=0){
long double offset=sqrtl(square);
if(clear(a.x+offset,side)||clear(a.x-offset,side)){
return true;
}
}
}
for(int j=0;j<i;j++){
Circle b=circles[j];
long double rb=b.radius+radius;
long double dx=b.x-a.x;
long double dy=b.y-a.y;
long double distance=hypotl(dx,dy);
if(distance==0||distance>ra+rb||distance<fabsl(ra-rb)){
continue;
}
long double along=(ra*ra-rb*rb+distance*distance)/(2*distance);
long double height=sqrtl(max(0.0L,ra*ra-along*along));
long double x=a.x+along*dx/distance;
long double y=a.y+along*dy/distance;
long double ox=-dy*height/distance;
long double oy=dx*height/distance;
if(clear(x+ox,y+oy)||clear(x-ox,y-oy)){
return true;
}
}
}
return false;
};
long double low=0;
long double high=1500;
for(int iteration=0;iteration<90;iteration++){
long double middle=(low+high)/2;
if(feasible(middle)){
low=middle;
}
else{
high=middle;
}
}
cout<<fixed<<setprecision(15)<<low<<'\n';
}