1012 Mow
这里我们可以分两类情况来讨论:人工除草便宜和机器除草便宜。如果人工除草便宜,答案直接为S*A(S为凸多边形面积,A为人工除草的单位面积价格)。
如果是机器除草便宜的话,先跑一个凸包,然后将凸包所有边向内缩小r个单位长度,求出新得到的凸多边形的核,并算出核的面积。总面积为:
S
核
=
S
核
+
C
核
∗
r
+
p
i
∗
r
2
S_核=S_核+C核*r+pi*r^2
S核=S核+C核∗r+pi∗r2
可以参考下图,其中蓝色多边形是原多边形(红色)的核:
#include<bits/stdc++.h>
#define db double
#define pb push_back
#define mem(a, b) memset(a, b, sizeof (a))
#define debug(x) cout << x << endl
using namespace std;
const db eps = 1e-10;
const db pi = acos(-1);
int sign(db x){ if(fabs(x) < eps) return 0; return x > 0 ? 1 : -1; }
int cmp(db k1, db k2){ return sign(k1-k2); }
int inmid(db k1, db k2, db k3){ return sign(k2-k1)*sign(k3-k1) <= 0; }//k1在k2、k3之间
struct point{
db x, y;
point(){}
point(db k1, db k2){ x = k1, y = k2; }
point operator + (const point &k1) const { return point(x+k1.x, y+k1.y); }//向量加法、点+向量=点
point operator - (const point &k1) const { return point(x-k1.x, y-k1.y); }//向量减法、点-点=向量
point operator * (db k1) const { return (point){x*k1, y*k1}; }//向量数乘
point operator / (db k1) const { return (point){x/k1, y/k1}; }//向量数除
int operator == (const point &k1) const { return cmp(x,k1.x)==0 && cmp(y,k1.y)==0; }//比较两个点(向量)是否相同
point turn(db k1){return (point){x*cos(k1)-y*sin(k1), x*sin(k1)+y*cos(k1)};}//逆时针旋转
point turn90(){return (point){-y, x};}//逆时针旋转90度
bool operator < (const point k1) const{
int a = cmp(x, k1.x);
if(a == -1) return 1; else if(a == 1) return 0; else return cmp(y,k1.y)==-1;
}//比较两个点(向量)的大小,x越小则点越小;若x相等,则y越小点越小。可以实现按点的坐标排序
db len(){ return sqrt(x*x+y*y); }//向量模长
db len2(){ return x*x+y*y; }//向量模长的平方
point unit(){ db k = len(); return point(x/k, y/k); }
db angle() { return atan2(y, x); }
point getdel(){if (sign(x)==-1||(sign(x)==0&&sign(y)==-1)) return (*this)*(-1); else return (*this);}
//当横坐标为负时,或横坐标为0纵坐标为负时,将点按原点做对称??? 角度是[-π/2, π/2]
int getp() const { return sign(y)==1 || (sign(y)==0&&sign(x)==-1); }
//判断点在12象限,或者在x的负半轴上??? 角度是(0, π]
void scan(){ scanf("%lf%lf", &x, &y); };
void print(){ printf("%.11f %.11f\n", x, y); }
};
int inmid(point k1, point k2, point k3){ return inmid(k1.x,k2.x,k3.x) && inmid(k1.y,k2.y,k3.y); }//k3 在 [k1,k2] 内
db dis2(point k1, point k2){ return (k1.x-k2.x)*(k1.x-k2.x) + (k1.y-k2.y)*(k1.y-k2.y); }
db dis(point k1, point k2){ return sqrt(dis2(k1, k2)); }
db cross(point k1, point k2){ return k1.x*k2.y - k1.y*k2.x; }//叉乘
db dot(point k1, point k2){ return k1.x*k2.x + k1.y*k2.y; }//点乘
db rad(point k1, point k2){ return atan2(cross(k1,k2), dot(k1,k2)); }//向量夹角
int compareangle (point k1,point k2){//极角排序,[-π, π]
return k1.getp()<k2.getp() || (k1.getp()==k2.getp() && sign(cross(k1,k2))>0);
}
point proj(point q, point k1, point k2){ point k=k2-k1; return k1+k*(dot(q-k1,k)/k.len2()); }// q 到直线 k1,k2 的投影
point reflect(point q, point k1, point k2){ return proj(q,k1,k2)*2-q; }// q 关于直线 k1, k2 的对称点
int clockwise(point k1,point k2,point k3){ return sign(cross(k2-k1,k3-k1)); }// k1 k2 k3 逆时针 1 顺时针 -1 否则 0
int checkll(point k1,point k2,point k3,point k4){// 直线相交判定 叉积计算面积,两直线不平行必相交(除去重合的情况),平行时,三角形面积相等
return cmp(cross(k3-k1,k4-k1),cross(k3-k2,k4-k2))!=0;
}
point getll(point k1,point k2,point k3,point k4){// 直线交点
db w1=cross(k1-k3,k4-k3),w2=cross(k4-k3,k2-k3); return (k1*w2+k2*w1)/(w1+w2);
}
int intersect(db l1,db r1,db l2,db r2){// 区间相交判定,区间1在区间2前
if(l1>r1) swap(l1,r1); if(l2>r2) swap(l2,r2); return cmp(r1,l2)!=-1 && cmp(r2,l1)!=-1;
}
db displ(point q, point k1, point k2){// 点到直线的距离
return fabs(cross(k2-k1, q-k1)) / (k2-k1).len();
}
struct line{// p[0]->p[1]
point p[2];
line(){}
line(point k1,point k2){p[0]=k1; p[1]=k2;}
point& operator [] (int k){return p[k];}
int include(point k){ return sign(cross(p[1]-p[0],k-p[0]))>0; }//点在直线左侧的判定,不包括在直线上
point dir(){return p[1]-p[0];}//方向向量
line push(db len) { // 向外(左手边)平移 len 个单位
point delta = (p[1] - p[0]).turn90().unit() * len;
return line(p[0] - delta, p[1] - delta);
}
};
int checkll(line k1,line k2){ return checkll(k1[0],k1[1],k2[0],k2[1]); }//直线相交判定
point getll(line k1,line k2){ return getll(k1[0],k1[1],k2[0],k2[1]); }//直线交点
bool checkpos(line k1, line k2, line k3) { return k3.include(getll(k1, k2)); }
int parallel(line k1,line k2){ return sign(cross(k1.dir(),k2.dir()))==0; }//直线平行判定,可以重合
bool sameDir(line k1, line k2) { return parallel(k1, k2) && sign(dot(k1.dir(), k2.dir())) == 1; }
bool operator < (line k1, line k2) {
if (sameDir(k1, k2)) return k2.include(k1[0]);
return compareangle(k1.dir(), k2.dir());
}
vector<point>poly;
int n;
point midpoint(point k1, point k2){//求中点
return point((k1.x+k2.x)/2, (k1.y+k2.y)/2);
}
line getaxis(point k1, point k2){//求k1k2的中垂线
point a1 = midpoint(k1, k2);
point a2 = a1 + (k1-k2).turn90();
return line(a1, a2);
}
bool online(point k1, line l){//点在直线上的判定
return sign(cross(l[1]-l[0], l[1]-k1)) == 0;
}
bool symmetry(point k1, point k2, line l){//k1k2关于l对称的判定
return sign(dot(k2-k1, l[1]-l[0])) == 0 && online(midpoint(k1, k2), l);
}
bool check(int x, int y){
line axis = getaxis(poly[x], poly[y]);
vector<point>pl, pr;
int change = 0, l, r;
l = x, r = y;
pl.pb(poly[x]);
pr.pb(poly[y]);
while(true){//把多边形上的点分成两部分
if(l == r){
if(!online(poly[l], axis)) change++;
break;
}
if(l != x && r != y) {
pl.pb(poly[l]);
pr.pb(poly[r]);
}
if((l-1+n)%n == r) break;
l = (l-1+n)%n, r = (r+1)%n;
}
l = (x+1)%n, r = (y-1+n)%n;
if(l == r && !online(poly[l], axis)) change++;//当check的点是i和i+2时,需要检查i+1这个点是否在对称轴上
int cnt = pl.size();
for(int i = 0; i < cnt; i++){//两部分中对应的点不能交叉
if(sign(cross(axis[1]-axis[0], pl[0]-axis[0])) != sign(cross(axis[1]-axis[0], pl[i]-axis[0])) &&
sign(cross(axis[1]-axis[0], pr[0]-axis[0])) != sign(cross(axis[1]-axis[0], pr[i]-axis[0])) ) return false;
}
for(int i = 0; i < cnt; i++){//判断需要移动的点的个数
if(!symmetry(pl[i], pr[i], axis)){
if(online(pl[i], axis) && online(pr[i], axis)) return false;
change++;
}
}
return change <= 1;
}
vector<point> convexHull(vector<point>A, int flag = 1) { // flag=0 不严格 flag=1 严格 //将坐标系存储变为逆时针
int n = A.size(); vector<point>ans(n + n);
sort(A.begin(), A.end()); int now = -1;
for (int i = 0; i < A.size(); i++) {
while (now > 0 && sign(cross(ans[now] - ans[now - 1], A[i] - ans[now - 1])) < flag) now--;
ans[++now] = A[i];
}
int pre = now;
for (int i = n - 2; i >= 0; i--) {
while (now > pre&& sign(cross(ans[now] - ans[now - 1], A[i] - ans[now - 1])) < flag)
now--;
ans[++now] = A[i];
}
ans.resize(now);
return ans;
}
db polygonArea(vector<point> P) {//求凸多边形面积
point top = P[0];
db ans = 0.0;
for (int i = 1; i < P.size() - 1; ++i)
ans += cross(P[i] - top, P[i + 1] - P[i]);
return fabs(0.5 * ans);
}
db polygonCircumference(vector<point> p) {//求凸多边形周长
db ans = 0.0;
int n = p.size();
for (int i = 0; i < n; ++i)
ans += fabs(dis(p[i],p[(i + 1) % n]));
return ans;
}
vector<line> getHL(vector<line>& L) { // 求半平面交 逆时针方向存储
sort(L.begin(), L.end());
deque<line> q;
for (int i = 0; i < (int)L.size(); ++i) {
if (i && sameDir(L[i], L[i - 1])) continue;
while (q.size() > 1 && !checkpos(q[q.size() - 2], q[q.size() - 1], L[i])) q.pop_back();
while (q.size() > 1 && !checkpos(q[1], q[0], L[i])) q.pop_front();
q.push_back(L[i]);
}
while (q.size() > 2 && !checkpos(q[q.size() - 2], q[q.size() - 1], q[0])) q.pop_back();
while (q.size() > 2 && !checkpos(q[1], q[0], q[q.size() - 1])) q.pop_front();
vector<line> ans;
for (int i = 0; i < q.size(); ++i) ans.push_back(q[i]);
return ans;
}
int main(){
ios::sync_with_stdio(false);
cin.tie(nullptr); cout.tie(nullptr);
int t; cin >> t;
while (t--) {
int n, r, a, b;
cin >> n >> r >> a >> b;
vector<point> p(n);
for (auto& i : p) cin >> i.x >> i.y;
p = convexHull(p);
//n = (int)p.size();
db s = polygonArea(p);
if (a <= b) {
db ans = s * a;
cout << fixed << setprecision(12) << ans << '\n';
}
else {
vector<line> l(n);
for (int i = 0; i < n; ++i) {
l[i] = line(p[i], p[(i + 1) % n]);
l[i] = l[i].push(-r);
}
vector<line> res = getHL(l);
if ((int)res.size() < 3) {
db ans = a * s;
cout << fixed << setprecision(12) << ans << '\n';
continue;
}
vector<point> pp((int)res.size());
for (int i = 0; i < (int)res.size(); ++i)
pp[i] = getll(res[i], res[(i + 1) % (int)res.size()]);
db ps = polygonArea(pp) + pi * r * r + polygonCircumference(pp) * r;
db ans = (s - ps) * a + ps * b;
cout << fixed << setprecision(12) << ans << '\n';
}
}
return 0;
}