[JZOJ5036]原谅

题目大意

平面上有$n$个点，每个点都有概率$p$出现（所有点概率一样）。
你可以在任意两个点之间连一条线段，两条线段不能在除了端点以外的地方相交。
你想要知道最多可以连接的线段的期望数目。
答案对$10^8+7$取模，$p$给定时已经取模。
数据保证没有三点共线。

$1\le n\le 10^3$

---

题目分析

显然对于任意一个点集，三角剖分一定是最优的答案。
根据欧拉公式，平面点集的任意一个三角剖分的边数都可以由下式决定：

$$3V-d-3=E$$
其中， $V$是点数，$d$是凸包上边数， $E$是三角剖分边数。
那么我们显然可以分开来求这条式子的期望。
$E[3V]=3E[V]=3np$。求 $E[d]$是经典问题：直接枚举一条边的起点，然后按照极角序枚举终点，使用指针来维护某个方向上的点数，限定某个方向不能有点来统计这条边出现在凸包上的概率。
注意，凸包只有一或两个点的情况（凸包边数分别为 $0$和$2$），按照上面的方法是不用特殊处理的。但是一个点都不选的话，我们根据上面的公式会算出 $-3$，最后乘上概率加回去就好了。
时间复杂度 $\mathrm O(n^2 \log n)$。

---

代码实现

#include <algorithm>
#include <iostream>
#include <cstdio>
#include <cctype>
#include <cmath>

using namespace std;

typedef long double db;

const int MOD=100000007;
const int N=1005;
const db EPS=1e-8;

bool equ(db x,db y){return fabs(x-y)<=EPS;}
int sgn(db x){return equ(x,0.0)?0:(x<0?-1:1);}
int sqr(int x){return 1ll*x*x%MOD;}

struct P
{
    db x,y;

    P (db x_=0.,db y_=0.){x=x_,y=y_;}

    P operator+(P const p)const{return P(x+p.x,y+p.y);}
    P operator-(P const p)const{return P(x-p.x,y-p.y);}
    P operator*(db const k)const{return P(x*k,y*k);}

    db operator*(P const p)const{return x*p.x+y*p.y;}
    db operator^(P const p)const{return x*p.y-y*p.x;}

    db arg(){return atan2(y,x);}
}pts[N],pt[N];

bool operator<(P p,P q){return p.arg()<q.arg();}

int pw[N],upw[N];
int n,p,up,Ev,Ed,ans,cnt;

void pre()
{
    pw[0]=1,upw[0]=1;
    for (int i=1;i<=n;++i) pw[i]=1ll*pw[i-1]*p%MOD,upw[i]=1ll*upw[i-1]*up%MOD;
    Ev=1ll*n*p%MOD;
}

void calc()
{
    for (int i=1;i<=n;++i)
    {
        cnt=0;
        for (int j=1;j<=n;++j)
            if (i!=j) pt[++cnt]=pts[j]-pts[i];
        sort(pt+1,pt+1+cnt);
        for (int cur=1,j=1,sign=0,tot;j<=cnt;++j)
        {
            if (!sign&&sgn(pt[j].arg())>=0) sign=1,cur=1;
            for (;cur<=cnt&&sgn(pt[j]^pt[cur])>=0&&(!sign||cur<j);++cur);
            tot=sign?cnt-j+cur-1:cur-j-1;
            (Ed+=1ll*upw[tot]*sqr(p)%MOD)%=MOD;
        }
    }
}

int main()
{
    freopen("forgive.in","r",stdin),freopen("forgive.out","w",stdout);
    scanf("%d%d",&n,&p),up=(1-p+MOD)%MOD;
    for (int i=1;i<=n;++i)
    {
        double x,y;
        scanf("%lf%lf",&x,&y),pts[i]=P(x,y);
    }
    pre(),calc(),ans=(((Ev*3-Ed-3)%MOD+MOD)%MOD+3*upw[n]%MOD)%MOD;
    printf("%d\n",ans);
    fclose(stdin),fclose(stdout);
    return 0;
}