NEUOJ 1210(The number of triangulation-计算几何+剪枝)

一个正n(<=400)边形分成n-3个三角形，面积最大的那个三角形只出现一次，求方案数。

如果没有面积最大的那个三角形只出现一次就是Catalan数
现在考虑枚举面积最大三角形，计算面积不超过S的正n边形的l条边组成的(l+1)边形的三角形划分方案数$dp_i$，暴力统计O(n^4)
加几个剪枝：
我们显然可以把全等的三角形一起计算。
但三角形较少的时候$dp_i=C_{i-1}$，$C_{i}$为Catalan数
预处理多边形所有顶点构成三角形的大小。
计算$dp_i=\sum_{j=2}^{i-1}dp_j*dp_{i-j+1}，只需要计算前一半$

另外这题是原题……

#include<bits/stdc++.h>
using namespace std;
#define For(i,n) for(int i=1;i<=n;i++)
#define Fork(i,k,n) for(int i=k;i<=n;i++)
#define Rep(i,n) for(int i=0;i<n;i++)
#define ForD(i,n) for(int i=n;i;i--)
#define ForkD(i,k,n) for(int i=n;i>=k;i--)
#define RepD(i,n) for(int i=n;i>=0;i--)
#define Forp(x) for(int p=Pre[x];p;p=Next[p])
#define Forpiter(x) for(int &p=iter[x];p;p=Next[p])  
#define Lson (o<<1)
#define Rson ((o<<1)+1)
#define MEM(a) memset(a,0,sizeof(a));
#define MEMI(a) memset(a,127,sizeof(a));
#define MEMi(a) memset(a,128,sizeof(a));
#define INF (2139062143)
#define pb push_back
#define mp make_pair 
#define fi first
#define se second
#define vi vector<int> 
#define pi pair<int,int>
#define SI(a) ((a).size())
#define Pr(kcase,ans) printf("Case #%d: %I64d\n",kcase,ans);
#define PRi(a,n) For(i,n-1) cout<<a[i]<<' '; cout<<a[n]<<endl;
#define PRi2D(a,n,m) For(i,n) { \
                        For(j,m-1) cout<<a[i][j]<<' ';\
                        cout<<a[i][m]<<endl; \
                        } 
#define ALL(x) (x).begin(),(x).end()
#pragma comment(linker, "/STACK:102400000,102400000")
typedef long long ll;
typedef long double ld;
typedef unsigned long long ull;
int read()
{
    int x=0,f=1; char ch=getchar();
    while(!isdigit(ch)) {if (ch=='-') f=-1; ch=getchar();}
    while(isdigit(ch)) { x=x*10+ch-'0'; ch=getchar();}
    return x*f;
} 
double sqr(double a){return a*a;}
const double eps=1e-10;
int dcmp(double x) {
    if (fabs(x)<eps) return 0; else return x<0 ? -1 : 1; 
}
ld PI = 3.141592653589793238462643383;
class P{
public:
    double x,y;
    P(double x=0,double y=0):x(x),y(y){}
    friend ld dis2(P A,P B){return sqr(A.x-B.x)+sqr(A.y-B.y);   }
    friend ld Dot(P A,P B) {return A.x*B.x+A.y*B.y; }
    friend ld Length(P A) {return sqrt(Dot(A,A)); }
    friend ld Angle(P A,P B) {
        if (dcmp(Dot(A,A))==0||dcmp(Dot(B,B))==0||dcmp(Dot(A-B,A-B))==0) return 0;
        return acos(max((ld)-1.0, min((ld)1.0, Dot(A,B) / Length(A) / Length(B) )) ); 
    }

    friend P operator- (P A,P B) { return P(A.x-B.x,A.y-B.y); }
    friend P operator+ (P A,P B) { return P(A.x+B.x,A.y+B.y); }
    friend P operator* (P A,double p) { return P(A.x*p,A.y*p); }
    friend P operator/ (P A,double p) { return P(A.x/p,A.y/p); }
    friend bool operator< (const P& a,const P& b) {return dcmp(a.x-b.x)<0 ||(dcmp(a.x-b.x)==0&& dcmp(a.y-b.y)<0 );}

}; 
P read_point() {
    P a;
    scanf("%lf%lf",&a.x,&a.y);
    return a;   
} 
bool operator==(const P& a,const P& b) {
    return dcmp(a.x-b.x)==0 && dcmp(a.y-b.y) == 0;
} 
typedef P V;

double Cross(V A,V B) {return A.x*B.y - A.y*B.x;}
double Area2(P A,P B,P C) {return fabs(Cross(B-A,C-A));}
V Rotate(V A,double rad) {
    return V(A.x*cos(rad)-A.y*sin(rad),A.x*sin(rad)+A.y*cos(rad));
} 
V Normal(V A) { 
    double L = Length(A);
    return V(-A.y/L , A.x/L); 
}

namespace complex_G{
    typedef complex<double> Point;
    typedef Point Vector;
    double Dot(Vector A,Vector B) {return real(conj(A)*B); }
    double Cross(Vector A,Vector B) {return imag(conj(A)*B); }
    Vector Rotate(Vector A,double rad) {return A*exp(Point(0,rad)); }
}
P GetLineIntersection(P p,V v,P Q,V w){
    V u = p-Q;
    double t = Cross(w,u)/Cross(v,w);
    return p+v*t;
}
P GetLineIntersectionB(P p,V v,P Q,V w){
    return GetLineIntersection(p,v-p,Q,w-Q);
}

double DistanceToLine(P p,P A,P B) {
    V v1 = B-A, v2 = p-A;
    return fabs(Cross(v1,v2))/Length(v1);
}
double DistanceToSegment(P p,P A,P B) {
    if (A==B) return Length(p-A);
    V v1 = B-A, v2 = p-A, v3 = p - B;
    if (dcmp(Dot(v1,v2))<0) return Length(v2);
    else if (dcmp(Dot(v1,v3))>0 ) return Length(v3);
    else return fabs(Cross(v1,v2) ) / Length(v1);
}
P GetLineProjection(P p,P A,P B) {
    V v=B-A;
    return A+v*(Dot(v,p-A)/Dot(v,v));
}
bool SegmentProperIntersection(P a1,P a2,P b1,P b2) { 
    double  c1 = Cross(a2-a1,b1-a1) , c2 = Cross(a2-a1,b2-a1),
            c3 = Cross(b2-b1,a1-b1) , c4 = Cross(b2-b1,a2-b1);
    return dcmp(c1)*dcmp(c2)<0 && dcmp(c3)*dcmp(c4)<0;
}
bool OnSegment(P p,P a1,P a2) {
    return dcmp(Cross(a1-p,a2-p)) == 0 && dcmp(Dot(a1-p,a2-p))<0;
}
double PolygonArea(P *p,int n) {
    double area=0;
    For(i,n-2) area+=Cross(p[i]-p[0],p[i+1]-p[0]);
    return area/2;
} 
void print(P p) {
    printf("(%.6lf,%.6lf)",p.x,p.y);
}
double torad(double deg) {
    return deg/180*acos(-1);
}
#define MAXN (3110)
int n;
ll x,y;
double S;
ll dp[MAXN];
P a[MAXN];
ll C[MAXN];
ll prework(int n,ll F) {
    C[1]=C[0]=1%F;
    Fork(i,2,n) {
        C[i]=0;
        Rep(j,i) C[i]+=C[j]*C[i-j-1]%F;
        C[i]%=F;
    }
}

ll h[5000000],hcase[5000000]={},tt=1;
int ha(int x) {
    if (hcase[x]!=tt) return -1;
    return h[x];
}
int hset(int x,int y) {
    h[x]=y;
    hcase[x]=tt;
}
int id(int x,int y) {
    return (x*(411))+y;
}
double sz[MAXN][MAXN];
int t[3];
ll calc(ll F) {
    For(i,n) For(j,n) sz[i][j]=Area2(a[1],a[i],a[j]);
    ll ans=0;
    prework(n,F);
    tt++;
    int tot=0;
    Rep(t0,n) Fork(t1,t0,n-t0+1) {
        int i=1,j=i+t0-1,k=j+t1-1;

        t[0]=t0,t[1]=t1,t[2]=i+n-k+1;
        if(t[0]>t[1]||t[1]>t[2]) continue;

        S=sz[j-i+1][k-i+1];

        ++tot;
        {
            dp[2]=1;
            dp[3]=dcmp(sz[2][3]-S)==-1;
            int ml=max(max(j-i,k-j),i+n-k)+1;
            Fork(i,4,ml) {
                dp[i]=0;
                if (sz[i][(i+1)/2]+eps<S) {
                    dp[i]=C[i-2]; continue;
                }
                Fork(j,2,i/2) 
                    if (dcmp(sz[i][j]-S)==-1)
                        (dp[i]+=2*dp[j]*dp[i-j+1]%F);
                    else break;
                if (i&1) {
                    int j=(i+1)/2;
                    if (dcmp(sz[i][j]-S)==-1)
                        (dp[i]+=dp[j]*dp[i-j+1]%F);
                }
                dp[i]%=F;
            }
        }

        ll c=dp[j-i+1]*dp[k-j+1]%F*dp[i+n-k+1]%F;
        if (t[0]==t[2]) (c*=n/3)%=F;
        else if (t[0]==t[1]||t[1]==t[2]) (c*=n)%=F;
        else (c*=n*2)%=F;
        (ans+=c)%=F;
    }
//  cout<<tot<<endl;
    return ans;
}
int main()
{
//  freopen("G.in","r",stdin);
//  freopen(".out","w",stdout);
    int T=read();
    while(T--) {
        n=read(); 
        x=read(),y=read();
        double ang=(n-2)/(double)n *acos(-1);
        a[0]=P(0,0); V v(1,0);
        For(i,n) a[i]=a[i-1]+v,v=Rotate(v,acos(-1)-ang);
        ll ans=calc(x)*calc(y);
        printf("%lld\n",ans);

    }
    return 0;
}