好像又鸽子了几天博客
poj真是神奇的网站……
基础知识
点积
(a·b=a.x*b.x+a.y*b.y)
$a$在$b$上的投影乘以$b$的模长
叉积
(a×b=a.x*b.y-a.y*b.x)
$a,b$围成的平行四边形的有向面积
判断线段相交（跨立实验）
由一条线段的端点向另一条线段两端点做叉积，判断符号是否一致
向量旋转
逆时针旋转$k$度：
inline vector turn (vector a,double k)
{
return vector(a.x*cos(k)-a.y*sin(k),a.x*sin(k)+a.y*cos(k));
}
判断点在多边形内部
由一点引一条射线，若与多边形有奇数个交点则在多边形内部，否则在多边形外
二维凸包
求凸包
二维平面上给定$n$个点求凸包
显然$y$最小的点应该在凸包上
我们把$y$最小的点先选出来，将其他点极角排序
按照极角序将点入栈，利用斜率单调性判断是否弹栈
#include<iostream>
#include<cstdio>
#include<algorithm>
#include<cmath>
using namespace std;
namespace red{
#define y1 qwq
#define int long long
#define eps (1e-10)
inline int read()
{
int x=0;char ch,f=1;
for(ch=getchar();(ch<'0'||ch>'9')&&ch!='-';ch=getchar());
if(ch=='-') f=0,ch=getchar();
while(ch>='0'&&ch<='9'){x=(x<<1)+(x<<3)+ch-'0';ch=getchar();}
return f?x:-x;
}
const int N=1e5+10;
int n;
int st[N],top;
double ret;
struct node
{
double x,y;
node (double tx=0,double ty=0){x=tx,y=ty;}
inline double operator ^ (const node &t) const
{
return x*t.y-y*t.x;
}
inline node operator - (const node &t) const
{
return node(x-t.x,y-t.y);
}
}a[N];
inline double sqr(double x){return x*x;}
inline double dis(node a,node b)
{
return sqrt(sqr(b.x-a.x)+sqr(b.y-a.y));
}
inline bool cmp1(node a,node b)
{
return a.x==b.x?a.y<b.y:a.x<b.x;
}
inline bool cmp2(node q,node b)
{
double tmp=(q-a[1])^(b-a[1]);
if(!tmp) return q.x==b.x?q.y>b.y:q.x<b.x;
return tmp>0;
}
inline void main()
{
n=read();
for(int i=1;i<=n;++i)
{
scanf("%lf%lf",&a[i].x,&a[i].y);
}
sort(a+1,a+n+1,cmp1);
int tot=0;
a[0].x=-1e9+7;
for(int i=1;i<=n;++i)
{
if(a[i].x!=a[tot].x||a[i].y!=a[tot].y) a[++tot]=a[i];
}
sort(a+2,a+tot+1,cmp2);
st[++top]=1;st[++top]=2;
for(int i=3;i<=tot;++i)
{
while(top>2&&((a[st[top-1]]-a[st[top]])^(a[i]-a[st[top]]))>=0) --top;
st[++top]=i;
}
for(int i=1;i<top;++i)
{
ret+=dis(a[st[i]],a[st[i+1]]);
}
ret+=dis(a[st[top]],a[st[1]]);
printf("%.2lf",ret);
}
}
signed main()
{
red::main();
return 0;
}
二维凸包面积
利用叉积
inline double cross(node a,node b,node c)
{
return (b-a)^(c-a);//叉积
}
inline double are()
{
double ret=0;
for(int i=1;i<n;++i) ret+=fabs(cross(a[1],a[i],a[i+1]));
return ret*2;
}
多次判断点在凸包内部
将凸包分成$n-2$个三角形，每次取中间的三角形，用叉积判断当前点在三角形的内部还是左侧或右侧，二分求解
没写过
动态二维凸包
初始有三个点，每次加入一个点，维护凸包
由于凸包只会扩大，所以以初始的三点围成的三角形的重心为基准点求其他点的极角，将点插入合适的位置，再判断两侧的点是否删除
用平衡树/$set$维护
注意每次插入一个点不一定只弹出相邻的两个点……
#include<iostream>
#include<cstdio>
#include<algorithm>
#include<cmath>
#include<set>
using namespace std;
namespace red{
#define y1 qwq
#define int long long
#define double long double
#define iter set<node>::iterator
#define eps (1e-8)
inline int read()
{
int x=0;char ch,f=1;
for(ch=getchar();(ch<'0'||ch>'9')&&ch!='-';ch=getchar());
if(ch=='-') f=0,ch=getchar();
while(ch>='0'&&ch<='9'){x=(x<<1)+(x<<3)+ch-'0';ch=getchar();}
return f?x:-x;
}
const int N=1e5+10;
int n,opt;
struct node
{
double x,y,ang;
node(double tx=0,double ty=0,double tang=0){x=tx,y=ty,ang=tang;}
friend bool operator < (node a,node b)
{
return (a.ang<b.ang)||(a.ang==b.ang&&a.x<b.x);
}
inline node operator + (const node b) const
{
return node(x+b.x,y+b.y,0);
}
inline double operator ^ (const node &b) const
{
return x*b.y-y*b.x;
}
inline node operator - (const node &b) const
{
return node(x-b.x,y-b.y,0);
}
}a[N],zx;
set<node> q;
iter it1,it2;
iter pre(iter x)
{
return x==q.begin()?--q.end():--x;
}
iter nxt(iter x)
{
return ++x==q.end()?q.begin():x;
}
inline void main()
{
n=read();
for(int i=1;i<=3;++i)
{
opt=read();
a[i].x=read(),a[i].y=read();
zx=zx+a[i];
a[i].x*=3,a[i].y*=3;
}
for(int i=1;i<=3;++i)
{
a[i].ang=atan2(a[i].y-zx.y,a[i].x-zx.x);
q.insert(a[i]);
}
for(int i=4;i<=n;++i)
{
opt=read();
a[i].x=read()*3,a[i].y=read()*3;
a[i].ang=atan2(a[i].y-zx.y,a[i].x-zx.x);
it2=q.lower_bound(a[i]);
if(it2==q.end()) it2=q.begin();
it1=pre(it2);
if(opt==1)
{
if(((*it2-a[i])^(*it1-a[i]))<=0) continue;
q.insert(a[i]);
iter now=q.find(a[i]);
iter tx=pre(now),ty=pre(tx);
while(((*ty-*now)^(*tx-*now))<=0) q.erase(tx),tx=ty,ty=pre(tx);
tx=nxt(now),ty=nxt(tx);
while(((*ty-*now)^(*tx-*now))>=0) q.erase(tx),tx=ty,ty=nxt(tx);
}
else
{
if(((*it2-a[i])^(*it1-a[i]))<=0) puts("YES");
else puts("NO");
}
}
}
}
signed main()
{
red::main();
return 0;
}
旋转卡壳
组合数学入门题
众所周知，旋转卡壳有$232*2=24$种读法（雾
其实就是单调性优化的枚举？
给定凸包求最大直径
然后我们发现枚举一个点的话另一个点有单调性，所以这个过程是$O(n)$的……
#include<iostream>
#include<cstdio>
#include<algorithm>
#include<cmath>
using namespace std;
namespace red{
#define y1 qwq
#define int long long
#define eps (1e-10)
inline int read()
{
int x=0;char ch,f=1;
for(ch=getchar();(ch<'0'||ch>'9')&&ch!='-';ch=getchar());
if(ch=='-') f=0,ch=getchar();
while(ch>='0'&&ch<='9'){x=(x<<1)+(x<<3)+ch-'0';ch=getchar();}
return f?x:-x;
}
const int N=1e5+10;
int n;
int st[N],top;
double ret;
struct node
{
double x,y;
node (double tx=0,double ty=0){x=tx,y=ty;}
inline double operator ^ (const node &t) const
{
return x*t.y-y*t.x;
}
inline node operator - (const node &t) const
{
return node(x-t.x,y-t.y);
}
}a[N];
inline double sqr(double x){return x*x;}
inline double dis(node a,node b)
{
return sqr(b.x-a.x)+sqr(b.y-a.y);
}
inline bool cmp1(node a,node b)
{
return a.x==b.x?a.y<b.y:a.x<b.x;
}
inline bool cmp2(node q,node b)
{
double tmp=(q-a[1])^(b-a[1]);
if(!tmp) return q.x==b.x?q.y>b.y:q.x<b.x;
return tmp>0;
}
inline double are(node a,node b,node c)
{
return fabs((a.x*b.y+b.x*c.y+c.x*a.y-a.x*c.y-b.x*a.y-c.x*b.y)/2);
}
inline void main()
{
n=read();
for(int i=1;i<=n;++i)
{
scanf("%lf%lf",&a[i].x,&a[i].y);
}
sort(a+1,a+n+1,cmp1);
int tot=0;
a[0].x=-1e9+7;
for(int i=1;i<=n;++i)
{
if(a[i].x!=a[tot].x||a[i].y!=a[tot].y) a[++tot]=a[i];
}
sort(a+2,a+tot+1,cmp2);
st[++top]=1;st[++top]=2;
for(int i=3;i<=tot;++i)
{
while(top>2&&((a[st[top-1]]-a[st[top]])^(a[i]-a[st[top]]))>=0) --top;
st[++top]=i;
}
tot=top;
for(int i=1;i<=tot;++i) st[++top]=st[i];
int tmp=3;
for(int i=1;i<=tot;++i)
{
while(are(a[st[i]],a[st[i+1]],a[st[tmp+1]])>are(a[st[i]],a[st[i+1]],a[st[tmp]])) ++tmp;
ret=max(ret,max(dis(a[st[i]],a[st[tmp]]),dis(a[st[i+1]],a[st[tmp]])));
}
printf("%lld",(int)ret);
}
}
signed main()
{
red::main();
return 0;
}
另一个经典应用：给定两凸包求两凸包上点的最近距离
取其中一个凸包的$y$最小的点，另一凸包$y$最大的点开始枚举，方向相反，这样就有单调性了
#include<cstdio>
#include<iostream>
#include<cmath>
#include<algorithm>
using namespace std;
namespace red{
#define int long long
#define ls(p) (p<<1)
#define rs(p) (p<<1|1)
#define eps (1e-9)
inline int read()
{
int x=0;char ch,f=1;
for(ch=getchar();(ch<'0'||ch>'9')&&ch!='-';ch=getchar());
if(ch=='-') f=0,ch=getchar();
while(ch>='0'&&ch<='9'){x=(x<<1)+(x<<3)+ch-'0';ch=getchar();}
return f?x:-x;
}
const int N=5e4+10;
int n,m;
struct node
{
double x,y;
node() {};
node (double _x,double _y):x(_x),y(_y){}
inline double operator ^ (const node &t) const {return (x*t.y-y*t.x);}
inline double operator * (const node &t) const {return (x*t.x+y*t.y);}
inline node operator - (const node &t) const {return node(x-t.x,y-t.y);}
inline bool operator < (const node &b) const {return atan2(y,x)<atan2(b.y,b.x);}
}a[N],b[N];
inline double sqr(double x) {return x*x;}
inline double dist(node a,node b)
{
return sqrt(sqr(a.x-b.x)+sqr(a.y-b.y));
}
double multi(node A,node B,node C)//点积
{
return (B-A)*(C-A);
}
double cross(node A,node B,node C)//叉积
{
return (B-A)^(C-A);
}
inline double disline(node A,node B,node C)
{
if(dist(A,B)<eps) return dist(B,C);
if(multi(A,B,C)<-eps) return dist(A,C);//用点积判断是否在线段内部
if(multi(B,A,C)<-eps) return dist(B,C);
return fabs(cross(A,B,C)/dist(A,B));//点到线段距离
}
double mindist(node A,node B,node C,node D)
{
return min(min(disline(A,B,C),disline(A,B,D)),min(disline(C,D,A),disline(C,D,B)));
}
inline void check(node a[],int n)
{
for(int i=1;i<n-1;++i)
{
double tmp=cross(a[i],a[i+1],a[i+2]);
if(tmp<-eps) return;
else if(tmp>eps)
{
reverse(a+1,a+n+1);
return;
}
}
}
inline double solve(node a[],node b[],int n,int m)
{
int ymina=1,ymaxb=1;
for(int i=1;i<=n;++i) if(a[i].y<a[ymina].y) ymina=i;
for(int i=1;i<=m;++i) if(a[i].y>a[ymaxb].y) ymaxb=i;
double tmp,ret=1e9+7;
a[n+1]=a[1],b[m+1]=b[1];
for(int i=1;i<=n;++i)
{
while(tmp=cross(a[ymina+1],b[ymaxb+1],a[ymina])-cross(a[ymina+1],b[ymaxb],a[ymina])<eps)
{
++ymaxb;
if(ymaxb>m) ymaxb=1;
}
ret=min(ret,mindist(a[ymina],a[ymina+1],b[ymaxb],b[ymaxb+1]));
++ymina;
if(ymina>n) ymina=1;
}
return ret;
}
inline void main()
{
while("haku")
{
n=read(),m=read();
if(!(n|m)) return ;
for(int i=1;i<=n;++i) scanf("%lf%lf",&a[i].x,&a[i].y);
for(int i=1;i<=m;++i) scanf("%lf%lf",&b[i].x,&b[i].y);
check(a,n);check(b,m);
printf("%.5lf ",min(solve(a,b,n,m),solve(b,a,m,n)));
}
}
}
signed main()
{
red::main();
return 0;
}
半平面交
有点恶心的东西
给定$n$条直线，每条直线保留左侧部分，求最后留下的面积
首先我们有$n^2$的暴力算法：每加入一条直线，就枚举以前的直线
如果在左侧则保留
如果在右侧把它扔掉
如果相交保留左侧端点，右侧改为交点
当然一般来说$n^2$是不够用的
我们可以先将所有直线按照与$x$轴的夹角排序，夹角相同的保留靠左的
用双端队列维护半平面交
维护方式：先保证队列中至少有两个元素
如果队尾的两个元素的交点在当前直线右侧，把队尾弹出
队首同理
但是记住一定要先弹出队尾，因为新加入的元素与一开始加入的元素围成的面积更小
如在这张图中，如果标号就是加入队列的顺序，那么当$3$加入的时候，先判断$1$就会踢$1$，先判断$2$就会踢$2$，但是新加入的元素和队首围成的面积会更小，所以应该先踢队尾
#include<bits/stdc++.h>
using namespace std;
namespace red{
#define int long long
#define ls(p) (p<<1)
#define rs(p) (p<<1|1)
#define eps (1e-8)
inline int read()
{
int x=0;char ch,f=1;
for(ch=getchar();(ch<'0'||ch>'9')&&ch!='-';ch=getchar());
if(ch=='-') f=0,ch=getchar();
while(ch>='0'&&ch<='9'){x=(x<<1)+(x<<3)+ch-'0';ch=getchar();}
return f?x:-x;
}
const int N=5000;
int n,m,tot,sum,head,tail;
struct node
{
double x,y;
inline double operator ^ (const node &t) const{return x*t.y-y*t.x;}
inline node operator - (const node &t) const{return (node){x-t.x,y-t.y};}
inline node operator + (const node &t) const{return (node){x+t.x,y+t.y};}
inline node operator * (const double &t) const{return (node){x*t,y*t};}
}a[N],c[N];
struct seg
{
node a,b;
double k;
seg(){}
seg(const node &aa,const node &bb):a(aa),b(bb){k=atan2(b.y,b.x);}
inline bool operator < (const seg &t) const
{
return k<t.k;
}
}q[N],que[N];
inline double cross(node a,node b,node c)
{
return (b-a)^(c-a);
}
node get_node(seg A,seg B)//求交点
{
node C=A.a-B.a;
double t=(B.b^C)/(A.b^B.b);
return A.a+A.b*t;
}
inline int dcmp(double x)
{
if(fabs(x)<eps) return 0;
return x>0?1:-1;
}
inline double are()
{
double are=0;
for(int i=head;i<tail;++i){ are+=fabs(cross(c[head],c[i],c[i+1])); }
return are/2;
}
inline void work()
{
head=tail=1;
que[tail]=q[1];
for(int i=2;i<=sum;++i)
{
while(head<tail&&(q[i].b^(c[tail-1]-q[i].a))<=eps) --tail;
while(head<tail&&(q[i].b^(c[head]-q[i].a))<=eps) ++head;
que[++tail]=q[i];
if(fabs(que[tail].b^que[tail-1].b)<=eps) //平行保留较左的
{
--tail;
if((que[tail].b^(q[i].a-que[tail].a))>eps) que[tail]=q[i];
}
if(head<tail) c[tail-1]=get_node(que[tail-1],que[tail]);//c[i]表示c[i]和c[i+1]的交点
}
while(head<tail&&(que[head].b^(c[tail-1]-que[head].a))<=eps) --tail;
if(tail-head<=1) return;
c[tail]=get_node(que[head],que[tail]);
}
inline void main()
{
n=read();
for(int i=1;i<=n;++i)
{
m=read();
for(int j=1;j<=m;++j)
{
a[j].x=read(),a[j].y=read();
}
for(int j=1;j<=m;++j)
{
++sum;
int t=j+1;
if(j==m) t=1;
q[sum]=seg(a[j],a[t]-a[j]);
}
}
sort(q+1,q+sum+1);
work();
printf("%.3lf",are());
}
}
signed main()
{
red::main();
return 0;
}
自适应辛普森积分法
$simple$积分法
求积分$int_^frac{cx+d}{ax+b}dx$
先来说什么是辛普森积分法
假设有函数$g(x)=Ax^2+Bx+C$
(int_{a}^{b}g(x)dx)
微积分基本定理：
(=int_{a}^{b}frac{A}{3}(b^3-a^3)+frac{B}{2}(b^2-a^2)+C(b-a))
大力化简
(=frac{A}{3}(b-a)(a^2+ab+b^2)+frac{B}{2}(b+a)(b-a)+C(b-a))
(=frac{b-a}{6}(2A(a^2+ab+b^2)+3B(b+a)+6C))
(=frac{b-a}{6}(2Aa^2+2Aab+2Ab^2+3Ba+3Bb+6C))
(=frac{b-a}{6}(Aa^2+Ba+C+Ab^2+Bb+C+4A(frac{a+b}{2})^2+4B(frac{a+b}{2})+4C))
(=frac{b-a}{6}(g(a)+g(b)+4g(frac{a+b}{2})))
得到辛普森积分公式
(int_{a}^{b}g(x)dx=frac{b-a}{6}(g(a)+g(b)+4g(frac{a+b}{2})))
不过这是二次函数的积分公式，和我们上面那个式子有啥关系？
我们知道，一个复杂的函数可以近似拟合为一个二次函数
如果$f(x)$的拟合二次函数为$g(x)$，那么
(int_{a}^{b}f(x)dx≈frac{b-a}{6}(f(a)+f(b)+4f(frac{a+b}{2})))
当然拟合肯定存在误差，不过当定义域越小的时候，这个误差也就越小
我们可以把原函数分段，对于每一段求拟合函数$g(x)$的值再加起来，只要分得的段足够小就可以得到正确答案
但是分的太小会导致答案不正确， 太多会降低程序效率
这个时候我们可以让程序自适应
二分递归答案的段数和区间长度，精度达到要求就停止递归返回答案
#include<iostream>
#include<cstdio>
#include<algorithm>
#include<cmath>
#include<set>
using namespace std;
namespace red{
#define y1 qwq
inline int read()
{
int x=0;char ch,f=1;
for(ch=getchar();(ch<'0'||ch>'9')&&ch!='-';ch=getchar());
if(ch=='-') f=0,ch=getchar();
while(ch>='0'&&ch<='9'){x=(x<<1)+(x<<3)+ch-'0';ch=getchar();}
return f?x:-x;
}
double a,b,c,d,l,r;
inline double f(double x)
{
return (c*x+d)/(a*x+b);
}
inline double simpson(double l,double r)
{
double mid=(l+r)/2;
return (f(l)+4*f(mid)+f(r))*(r-l)/6;
}
inline double asr(double l,double r,double eps,double ans)
{
double mid=(l+r)/2;
double tl=simpson(l,mid),tr=simpson(mid,r);
if(fabs(tl+tr-ans)<=15*eps) return tl+tr+(tl+tr-ans)/15;
return asr(l,mid,eps/2,tl)+asr(mid,r,eps/2,tr);
}
inline double asr(double l,double r,double eps)
{
return asr(l,r,eps,simpson(l,r));
}
inline void main()
{
scanf("%lf%lf%lf%lf%lf%lf",&a,&b,&c,&d,&l,&r);
printf("%.6f
",asr(l,r,1e-6));
}
}
signed main()
{
red::main();
return 0;
}