[spoj CIRU2]關于圓面積并的思考，非simpson?。?/h1>

2019-11-10 17:04:12

字體：大 中 小

來源：轉載

供稿：網友

aekdycoin神犇的講解 圖一是普通求圓并的思路，復雜度O（n^2 logn） 圖二是神犇給的提示，沒有題解…

一開始困擾了好久，我覺得這和hzwer開的壞頭有關，他在BZOJ上圓并用了simpson之后一大批同學思路都跟著“奇怪”了起來，這題simpson沒試過，不過評論說會掛精度。 首先我們要理解非simpson的圓并在干什么？這個算法很神的地方在于他把圓并分成了好多個圓弧和一個多邊形，然后分別求面積，圓弧部分是一個圓上未被其它圓覆蓋的部分的圓?。赡苡卸喽危?，而多邊形則被分解為每條邊單獨算貢獻。這樣我們的工作就簡單許多了，枚舉一個圓，與其它圓求交點，搞出區間后排序，統計未覆蓋的角度算圓弧，統計坐標算叉積求貢獻。最后統統加起來就好了。 注意一些被包含的圓沒有作用，標記掉即可！ 現在再看那張提示圖，我們發現一個圓未被其它圓覆蓋的部分就是圓并的“一次覆蓋面積”（見spoj ciru2），我們大膽的猜一個結論：“k次覆蓋面積”就是k次覆蓋的圓弧加多邊形面積。呵呵...結果就是這樣... 注意這時被包含的圓不是去掉，是統計時覆蓋次數加一！ 廢話就這么多，總之細節很煩，對拍造數據要考慮數據的合法性...

普通圓并

#include <cmath>#include <cstdio>#include <cstdlib>#include <cstring>#include <iostream>#include <algorithm>using namespace std;typedef long long LL;inline int read() { int x=0,f=1; char ch=getchar(); while (ch<'0'||ch>'9') { if (ch=='-') f=-1; ch=getchar(); } while (ch>='0'&&ch<='9') { x=x*10+ch-'0'; ch=getchar(); } return x*f;}const int MAXN=1003;const double PI=3.141592653589793238462,eps=1e-8;int n,m;double ans1,ans2,b[MAXN],R[MAXN];struct P{ double x,y; P(double v=0,double u=0){ x=v; y=u; } double len(){ return sqrt(x*x+y*y); } P Operator +(const P &E)const{ return P(x+E.x,y+E.y); } P operator -(const P &E)const{ return P(x-E.x,y-E.y); } P operator *(const double &z)const{ return P(x*z,y*z); } P operator /(const double &z)const{ return P(x/z,y/z); } double operator *(const P &E)const{ return x*E.x+y*E.y; }}a[MAXN];inline double cross(P p,P q){ return p.x*q.y-p.y*q.x; }inline double dist(P p,P q){ return sqrt((p.x-q.x)*(p.x-q.x)+(p.y-q.y)*(p.y-q.y)); }struct Data{ P ax,ay; double x,y; bool operator <(const Data &E)const{ return (x<E.x)||((x==E.x)&&(y<E.y)); }}c[MAXN*2];inline P rot(P p,double a){ return P(cos(a)*p.x-sin(a)*p.y,sin(a)*p.x+cos(a)*p.y); }inline double deg(P p){ double tmp=atan2(p.y,p.x); return tmp<0?tmp+PI*2:tmp;}inline void work(P o1,P o2,double r1,double r2){ double s,p,h,alpha,tmp=dist(o1,o2); p=(r1+r2+tmp)/2; s=sqrt(p*(p-r1)*(p-r2)*(p-tmp)); h=s*2/tmp; alpha=asin(h/r1); if (r1*r1+tmp*tmp<r2*r2) alpha=PI-alpha; m++; c[m].ax=rot(o2-o1,PI*2-alpha)/tmp*r1+o1; c[m].ay=rot(o2-o1,alpha)/tmp*r1+o1; c[m].x=deg(c[m].ax-o1); c[m].y=deg(c[m].ay-o1); if (c[m].x>c[m].y) { m++; c[m].y=c[m-1].y; c[m].ay=c[m-1].ay; c[m-1].y=PI*2; c[m].x=0; c[m-1].ay=c[m].ax=o1+P(r1,0); }}int main() { //freopen("1.in","r",stdin); //freopen("2.out","w",stdout); int i,j,k,l; m=read(); for (i=1,n=0;i<=m;++i) { n++; a[n].x=read(); a[n].y=read(); R[n]=read(); if (fabs(R[n])<=eps) n--; } memset(b,0,sizeof(b)); for (i=1;i<=n;++i) if (!b[i]) for (j=1;j<=n;++j) if (!b[j]&&i!=j&&dist(a[i],a[j])<=R[j]-R[i]) { b[i]=1; break; } //for (i=1;i<=n;++i) if (!b[i]) spoj ciru2

#include <cmath>#include <cstdio>#include <cstdlib>#include <cstring>#include <iostream>#include <algorithm>using namespace std;typedef long long LL;inline int read() { int x=0,f=1; char ch=getchar(); while (ch<'0'||ch>'9') { if (ch=='-') f=-1; ch=getchar(); } while (ch>='0'&&ch<='9') { x=x*10+ch-'0'; ch=getchar(); } return x*f;}const int MAXN=1003;const double PI=3.1415926535897932384626433832795028841971693993,eps=1e-8;int n,m;double ans[MAXN][2],b[MAXN],R[MAXN];struct P{ double x,y; P(double v=0,double u=0){ x=v; y=u; } double len(){ return sqrt(x*x+y*y); } P operator +(const P &E)const{ return P(x+E.x,y+E.y); } P operator -(const P &E)const{ return P(x-E.x,y-E.y); } P operator *(const double &z)const{ return P(x*z,y*z); } P operator /(const double &z)const{ return P(x/z,y/z); } double operator *(const P &E)const{ return x*E.x+y*E.y; }}a[MAXN];inline double cross(P p,P q){ return p.x*q.y-p.y*q.x; }inline double dist(P p,P q){ return sqrt((p.x-q.x)*(p.x-q.x)+(p.y-q.y)*(p.y-q.y)); }struct Data{ P ax; double x; int y; bool operator <(const Data &E)const{ return x<E.x; }}c[MAXN*10];inline P rot(P p,double a){ return P(cos(a)*p.x-sin(a)*p.y,sin(a)*p.x+cos(a)*p.y); }inline double deg(P p){ double tmp=atan2(p.y,p.x); return tmp<0?tmp+PI*2:tmp;}inline void work(P o1,P o2,double r1,double r2){ double s,p,h,alpha,tmp=dist(o1,o2); p=(r1+r2+tmp)/2; s=sqrt(p*(p-r1)*(p-r2)*(p-tmp)); h=s*2/tmp; alpha=asin(h/r1); if (r1*r1+tmp*tmp<r2*r2) alpha=PI-alpha; c[++m].ax=rot(o2-o1,PI*2-alpha)/tmp*r1+o1; c[m].x=deg(c[m].ax-o1); c[m].y=1; c[++m].ax=rot(o2-o1,alpha)/tmp*r1+o1; c[m].x=deg(c[m].ax-o1); c[m].y=-1; if (c[m-1].x>c[m].x) { c[++m].ax=o1+P(r1,0); c[m].x=0; c[m].y=1; c[++m].ax=o1+P(r1,0); c[m].x=PI*2; c[m].y=-1; }}int main() { //freopen("1.in","r",stdin); //freopen("2.out","w",stdout); int i,j,k,l; n=read(); for (i=1;i<=n;++i) { a[i].x=read(); a[i].y=read(); R[i]=read(); } memset(b,0,sizeof(b)); for (i=1;i<=n;++i) for (j=1;j<=n;++j) if (i!=j&&dist(a[i],a[j])<=R[j]-R[i]) b[i]++; memset(ans,0,sizeof(ans)); for (i=1;i<=n;++i) { for (j=1,m=0;j<=n;++j) if (i!=j&&dist(a[i],a[j])<R[i]+R[j]&&dist(a[i],a[j])>fabs(R[j]-R[i])) work(a[i],a[j],R[i],R[j]); c[++m].ax=a[i]+P(R[i],0); c[m].x=0; c[m].y=b[i]+1; c[++m].ax=a[i]+P(R[i],0); c[m].x=PI*2; c[m].y=-b[i]-1; sort(c+1,c+1+m); /*printf("%d/n",i); for (j=1;j<=m;++j) { printf("%d %.3lf %d/n",j,c[j].x,c[j].y); printf("(%.3lf,%.3lf)/n",c[j].ax.x,c[j].ax.y); }*/ int sa=1,sb=0,sc; for (j=1;j<=m;j=k+1) { for (k=j;k<m&&fabs(c[k+1].x-c[j].x)<eps;++k); for (l=j,sc=sb;l<=k;++l) sc+=c[l].y; if (sc==sb&&k!=m) continue; if (fabs(c[sa].x-c[j].x)>eps) { ans[sb][0]+=cross(c[sa].ax,c[j].ax); double tmp=c[j].x-c[sa].x; ans[sb][1]+=R[i]*R[i]*(tmp-sin(tmp))/2; } sb=sc; sa=j; } } //for (i=1;i<=n;++i) printf("%.3lf %.3lf %.3lf/n",ans[i][0],ans[i][1],fabs(ans[i][0])/2+ans[i][1]); for (i=1;i<=n;++i) ans[i][0]=fabs(ans[i][0])/2+ans[i][1]; for (i=1;i<=n;++i) printf("[%d] = %.3lf/n",i,ans[i][0]-ans[i+1][0]); //fclose(stdin); //fclose(stdout); return 0;}/*53 0 34 3 51 2 44 3 51 4 4*/ 身為一個蒟蒻非常同情蒟蒻寫完程序死WA不止，希望我的程序對蒟蒻們的學習(對拍)有幫助。

上一篇：深入理解JVM(一)——JVM內存模型

下一篇：最大子段和

發表評論 共有條評論 用戶名: 密碼: 驗證碼: 匿名發表

學習交流

更多

• 

  索泰發布一款GTX 1070 Mini迷

• 

  AMD新旗艦顯卡輕松干翻NVIDIA

索泰發布一款GTX 1070 Mini迷你版本:小機

索泰發布一款GTX 1070 Mini迷你版本:小機箱大愛...

• usb無線網卡怎么用,小編告訴你安裝教程09-10
• usb調試在哪,小編告訴你usb調試在哪09-10
• 優盤不顯示,小編告訴你優盤不顯示怎么辦09-10
• 低級格式化,小編告訴你硬盤怎么低級格式化09-10

• 帝國cms分類信息的所在地在的修改09-08
• 將網站地圖和友情鏈接table樣式改為div+css09-08
• 用帝國cms實現不規則新聞或信息調用（應大站09-08
• 帝國調用DZ論壇精華帖09-08
• 用靈動標簽調用discuz和phpwind的最新貼子09-08

熱門圖片

更多

• 

  芭蕾舞蹈表演，真實美到極致

• 

  下午茶時間，悠然自得的休憩

• 

  充斥這繁華奢靡氣息的城市迪拜風景圖片

• 

  從山間到田野再到大海美麗的自然風景圖片

• 

  肉食主義者的最愛美食烤肉圖片

• 

  夏日甜心草莓美食圖片

• 

  人逢知己千杯少，喝酒搞笑圖集

• 

  搞笑試卷，學生惡搞答題

猜你喜歡的新聞

• 榮耀總裁趙明烏鎮演講：榮耀首款5G手機V30下
• 搜狐張朝陽：回歸媒體是搜狐重新崛起的關鍵
• 華為輪值董事長郭平：虛擬技術創造現實價值
• 第六屆世界互聯網大會開幕“to B”端成熱門
• 滴滴英文服務上線兩周年 用戶已超200萬
• 華為推出全球至快AI訓練集群Atlas900
• 馬斯克：特斯拉正組建中國技術團隊
• 10年后6G將問世 速度有望比5G快100倍
• WeworkCEO稱已開始考慮未來職位 不排除放棄
• 谷歌軟件商店模式變革：推出5美元會員 可用數

猜你喜歡的關注

• 【POJ 2528】Mayor&#39;s posters
• spring maven 搭建dubbo框架(dubbo-admin)
• 【POJ 3667】Hotel
• 【POJ 2104】K-th Number&amp;主席樹詳解
• flex4 spark 布局
• CUDA線程協作之共享存儲器“__shared__”&a
• 1038. Recover the Smallest Number (30)
• 字幕文件批量重命名腳本 —— Linux
• 1043. 輸出PATest(20)
• PAT甲級1003