[CF 559D] Randomizer

作法：

題目要求的是任取一個子多邊形的內部格點數期望值，顯然必須先使用 Pick's 定理，因為我們有$A=i+\frac{b}{2}-1$，所以我們只要知道子多邊形的$A$值期望值，和$b$值期望值就可以了。首先看$A$值的期望值，對於每個子多邊形，如果我們用整個多邊形去扣掉他，那麼他就會剩下一堆「頂點是原多邊形中的連續幾個點」的子多邊形。因此我們只要先算好要扣掉部份的期望值，再用原本的整個面積去扣掉就可以了。至於要如何算扣掉的部份，對於一個「頂點是原多邊形上連續幾個點」的子多邊形來說，假設他總共有$i$個頂點，那麼不難算出他總共會被算到$2^{n-i}-1$次，並且總共有$2^n-1-n-C_2^n$種子多邊形，因此有$i$個頂點的子多邊形貢獻為$\displaystyle \frac{2^{n-i}-1}{2^n-1-n-C_2^n}\times $其面積。這樣我們可以先獲得一個$O(n^2)$的算法：直接把每一項都加起來，並且在算「有$i$個頂點的子多邊形的面積總和」時可以用$i-1$個頂點的總和花$O(n)$的時間推過來。但觀察這條式子可以發現，當$i$越來越大時其倍率是指數往下跑的，因此我們可以只計算$i$很小的前幾項就好了，稍微估一下大概可以知道一直算到$i=100$左右就很夠了。

再來則是計算$b$的值，我們可以把每條邊拆開看，每條邊的貢獻會是$gcd(x$座標差$,y$座標差$)$（這樣算才能保證把子多邊形的每個邊的貢獻加起來是我們要的值），並且我們想知道他會被算到幾次。假設以這條邊切開原多邊形，所得到的子多邊形其中一邊有$k$個點，那麼不難得出他被算到的次數為：$(2^{n-2-k}-1)+(2^k-1)$，因此我們可以把他拆成兩邊看，也就是我們把主角放在「原多邊形的一段頂點區間」而不是「邊」，這樣就可以套用剛才的作法了，只要考慮那些「頂點區間夠小」的區間就可以了，也是取到長度$100$左右就很夠了。

code :

#include<bits/stdc++.h>
#define DB double
using namespace std;
const DB eps=1e-9 ;
const int maxn=100000+10 ;
 
int dcmp(DB x)
{
    if(fabs(x)<eps) return 0 ;
    return x>0 ? 1 : -1 ;
}
 
struct pt{DB x,y ;};
 
pt operator + (const pt &a,const pt &b) { return (pt){a.x+b.x,a.y+b.y} ; }
pt operator - (const pt &a,const pt &b) { return (pt){a.x-b.x,a.y-b.y} ; }
pt operator * (const pt &a,const DB &d) { return (pt){d*a.x,d*a.y} ; }
bool operator == (const pt &a,const pt &b) { return !dcmp(a.x-b.x) && !dcmp(a.y-b.y) ; }
bool operator != (const pt &a,const pt &b) { return !(a==b) ; }
 
DB dot(const pt &a,const pt &b) { return a.x*b.x + a.y*b.y ; }
DB cross(const pt &a,const pt &b) { return a.x*b.y - a.y*b.x ; }
DB cross(const pt &a,const pt &b,const pt &c) { return cross(b-a,c-a) ; }
DB dis2(const pt &a,const pt &b) { return dot(a-b,a-b) ; }
DB length(const pt &a) { return sqrt(dot(a,a)) ; }
 
typedef vector<pt> poly ;
DB Area(const poly &p)
{
    DB ret=0 ;
    for(int i=1;i+1<p.size();i++)
        ret+=cross(p[0],p[i],p[i+1])/2 ;
    return fabs(ret) ;
}
inline DB Area(const pt &a,const pt &b,const pt &c)
{
    return fabs(cross(a,b,c))/2 ;
}
 
int n,x[maxn],y[maxn] ;
poly p ;
DB pw2[maxn] ;
DB cal_area()
{
    DB now=0,sub=0 ;
    for(int i=1;i<=min(n-3,50);i++)
    {
        for(int j=0;j<n;j++) now+=Area(p[j],p[(j+i)%n],p[(j+i+1)%n]) ;
        sub+=now*(pw2[i+2]-pw2[n]) ;
    }
    sub*=(1+((n*n+n+2)/(pow(2,n+1)-(n*n+n+2)))) ;
    return Area(p)-sub ;
}
 
inline int getnum(int i,int j)
{
    int g=__gcd(x[i]-x[j],y[i]-y[j]) ;
    return g<0 ? -g : g ;
}
DB cal_edge()
{
    DB ret=0 ;
    for(int i=1;i<=min(n-2,50);i++)
        for(int j=0;j<n;j++)
        ret+=getnum(j,(j+i)%n)*(pw2[i+1]-pw2[n]) ;
    ret*=(1+((n*n+n+2)/(pow(2,n+1)-(n*n+n+2)))) ;
    return ret ;
}
 
main()
{
    scanf("%d",&n) ; p.resize(n) ;
    for(int i=0;i<maxn;i++) pw2[i]=(i ? pw2[i-1]*0.5 : 1) ;
    for(int i=0;i<n;i++)
    {
        scanf("%d%d",&x[i],&y[i]) ;
        p[i]=(pt){x[i],y[i]} ;
    }
    printf("%.15f\n",cal_area()-cal_edge()/2+1) ;
}