计算几何整理总结 以前学习的做的计算几何部分笔记，现在整理一下发出来 二维几何 预备部分 第一部分—点 第二部分—线 第三部分—圆 第四部分—多边形 其他 最小矩形面积覆盖 直线切凸多边形 半平面交 三维几何 平面最近点对 三维凸包 Sorted & Add By: MFDY From: kuangbin 二维几何 预备部分 const double eps = 1e-8; const double inf = 1e20; const double pi = acos(-1.0); const int maxp = 1010; //`Compares a double to zero` int sgn(double x) { if(fabs(x) < eps) return 0; if(x < 0) return -1; else return 1; } //square of a double inline double sqr(double x){return x*x;} 第一部分—点 point() - Empty constructor point(double _x,double _y) - constructor input() - double input output() - %.2f output operator == - 比较x和y operator < - 先比x，后比y operator - - x，y分别减; 从A点指向B点的矢量AB可用B-A来表示 operator ^ - 叉积:求两个矢量形成的面积;求面积，求顺逆时针方向，判断是否在半平面上 因为a*b为有向面积，可正可负，a逆时针旋转小于180为正 operator * - 矢量的点积,求投影长度 len() - gives length from origin len2() - gives square of length from origin distance(point p) - gives distance from p operator + point b - returns new point after adding curresponging x and y operator * double k - returns new point after multiplieing x and y by k operator / double k - returns new point after divideing x and y by k rad(point a,point b)- 计算pa和pb的夹角 trunc(double r) - 更改为长度为r的向量 rotleft() - 逆时针旋转90° rotright() - 顺时针旋转90° rotate(point p,double angle) - 绕p点逆时针旋转angle° struct Point{ double x, y; Point(){} Point(double _x, double _y){ x = _x; y = _y; } void input(){ scanf("%lf%lf", &x, &y); } void output(){ printf("%.2f %.2f\n", x, y); } bool operator == (Point b) const { return sgn(x - b.x) == 0 && sgn(y - b.y) == 0; } bool operator < (Point b) const { return sgn(x-b.x) == 0 ? sgn(y-b.y) < 0 : x < b.x; } Point operator -(const Point &b)const{ return Point(x-b.x,y-b.y); } //叉积 double operator ^(const Point &b)const{ return x*b.y - y*b.x; } //点积 double operator *(const Point &b)const{ return x*b.x + y*b.y; } //返回长度 double len(){ return hypot(x,y);//库函数 } //返回长度的平方 double len2(){ return x*x + y*y; } //返回两点的距离 double distance(Point p){ return hypot(x-p.x,y-p.y); } Point operator +(const Point &b)const{ return Point(x+b.x,y+b.y); } Point operator *(const double &k)const{ return Point(x*k,y*k); } Point operator /(const double &k)const{ return Point(x/k,y/k); } //`计算pa 和 pb 的夹角` //`就是求这个点看a,b 所成的夹角` //`测试 LightOJ1203` double rad(Point a,Point b){ Point p = *this; return fabs(atan2( fabs((a-p)^(b-p)),(a-p)*(b-p) )); } //`化为长度为r的向量` Point trunc(double r){ double l = len(); if(!sgn(l))return *this; r /= l; return Point(x*r,y*r); } //`逆时针旋转90度` Point rotleft(){ return Point(-y,x); } //`顺时针旋转90度` Point rotright() { return Point(y,-x); } //`绕着p点逆时针旋转angle` Point rotate(Point p,double angle){ Point v = (*this) - p; double c = cos(angle), s = sin(angle); return Point(p.x + v.x*c - v.y*s,p.y + v.x*s + v.y*c); } }; 第二部分—线 Stores two points line() - Empty constructor line(point _s,point _e) - line through _s and _e operator == - checks if two points are same line(point p,double angle) - 根据一个点和倾斜角angle确定直线，0≤angle<pi line(double a,double b,double c) - ax + by + c = 0 转化 input() - inputs s and e adjust() - orders in such a way that s < e length() - 求线段长度 angle() - 返回直线倾斜角 relation(point p) - 点和直线关系(向量) 1 点在直线左侧 2 点在直线右侧 3 点在直线上 xmult(Point p) - 点和直线叉积, 用于判断点和直线关系(优先选择！！！) pointonseg(double p) - 判断点是否在线段上 parallel(line v) - 判断两向量是否平行 segcrossseg(line v) - 两线段相交判断 returns 0 不相交 returns 1 非规范相交 returns 2 规范相交 linecrossseg(line v) - 直线和线段相交判断，同上 linecrossline(line v) - 两直线关系 0 平行 1 重合 2 相交 crosspoint(line v) -求两直线的交点,要保证两直线不平行或重合 dispointtoline(point p) - 点到直线的距离 dispointtoseg(point p) - 点到线段的距离 dissegtoseg(line v) - 返回线段到线段的距离,前提是两线段不相交，相交距离就是 0 了 lineprog(point p) - 返回点 p 在直线上的投影 symmetrypoint(point p) - 返回点 p 关于直线的对称点 struct Line{ Point s,e; Line(){} Line(Point _s,Point _e){ s = _s; e = _e; } bool operator ==(Line v){ return (s == v.s)&&(e == v.e); } //`根据一个点和倾斜角angle确定直线,0<=angle<pi` Line(Point p,double angle){ s = p; if(sgn(angle-pi/2) == 0){ e = (s + Point(0,1)); } else{ e = (s + Point(1,tan(angle))); } } //ax+by+c=0 Line(double a,double b,double c){ if(sgn(a) == 0){ s = Point(0,-c/b); e = Point(1,-c/b); } else if(sgn(b) == 0){ s = Point(-c/a,0); e = Point(-c/a,1); } else{ s = Point(0,-c/b); e = Point(1,(-c-a)/b); } } void input(){ s.input(); e.input(); } void adjust(){ if(e < s)swap(s,e); } //求线段长度 double length(){ return s.distance(e); } //`返回直线倾斜角 0<=angle<pi` double angle(){ double k = atan2(e.y-s.y,e.x-s.x); if(sgn(k) < 0)k += pi; if(sgn(k-pi) == 0)k -= pi; return k; } //点和直线关系(向量) //1 在左侧 //2 在右侧 //3 在直线上 int relation(Point p){ int c = sgn((p-s)^(e-s)); if(c < 0)return 1; else if(c > 0)return 2; else return 3; } // 点和直线的叉积，也可用于判断点和直线关系(优选) // return < 0 点在直线左侧, 即点在直线的逆时针方向 // return > 0 点在直线右侧, 即点在直线的顺时针方向 // return == 0 点在直线上 int xmult(Point p) { return (s - p) ^ (e - p); } // 点在线段上的判断 bool pointonseg(Point p){ return sgn((p-s)^(e-s)) == 0 && sgn((p-s)*(p-e)) <= 0; } //`两向量平行(对应直线平行或重合)` bool parallel(Line v){ return sgn((e-s)^(v.e-v.s)) == 0; } //`两线段相交判断` //`2 规范相交` //`1 非规范相交` //`0 不相交` int segcrossseg(Line v){ int d1 = sgn((e-s)^(v.s-s)); int d2 = sgn((e-s)^(v.e-s)); int d3 = sgn((v.e-v.s)^(s-v.s)); int d4 = sgn((v.e-v.s)^(e-v.s)); if( (d1^d2)==-2 && (d3^d4)==-2 )return 2; return (d1==0 && sgn((v.s-s)*(v.s-e))<=0) || (d2==0 && sgn((v.e-s)*(v.e-e))<=0) || (d3==0 && sgn((s-v.s)*(s-v.e))<=0) || (d4==0 && sgn((e-v.s)*(e-v.e))<=0); } //`直线和线段相交判断` //`-*this line -v seg` //`2 规范相交` //`1 非规范相交` //`0 不相交` int linecrossseg(Line v){ int d1 = sgn((e-s)^(v.s-s)); int d2 = sgn((e-s)^(v.e-s)); if((d1^d2)==-2) return 2; return (d1==0||d2==0); } //`两直线关系` //`0 平行` //`1 重合` //`2 相交` int linecrossline(Line v){ if((*this).parallel(v)) return v.relation(s)==3; return 2; } //`求两直线的交点` //`要保证两直线不平行或重合` Point crosspoint(Line v){ double a1 = (v.e-v.s)^(s-v.s); double a2 = (v.e-v.s)^(e-v.s); return Point((s.x*a2-e.x*a1)/(a2-a1),(s.y*a2-e.y*a1)/(a2-a1)); } //点到直线的距离 double dispointtoline(Point p){ return fabs((p-s)^(e-s))/length(); } //点到线段的距离 double dispointtoseg(Point p){ if(sgn((p-s)*(e-s))<0 || sgn((p-e)*(s-e))<0) return min(p.distance(s),p.distance(e)); return dispointtoline(p); } //`返回线段到线段的距离` //`前提是两线段不相交，相交距离就是0了` double dissegtoseg(Line v){ return min(min(dispointtoseg(v.s),dispointtoseg(v.e)),min(v.dispointtoseg(s),v.dispointtoseg(e))); } //`返回点p在直线上的投影` Point lineprog(Point p){ return s + ( ((e-s)*((e-s)*(p-s)))/((e-s).len2()) ); } //`返回点p关于直线的对称点` Point symmetrypoint(Point p){ Point q = lineprog(p); return Point(2*q.x-p.x,2*q.y-p.y); } }; 第三部分—圆 circle(Point a,Point b,Point c) - 三角形的外接圆 circle(Point a,Point b,Point c,bool t) - 三角形的内切圆 input() - 输入 output() - 输出 area() - 面积 circumference() - 周长 relation() - 点和圆的关系 圆外 圆上 圆内 relationseg(Line v) - 线段和圆的关系 relationline(Line v) - 直线和圆的关系 relationcircle(circle v) - 两圆的关系 pointcrosscircle(circle v,Point &p1,Point &p2) - 两个圆的交点个数 pointcrossline(Line v,Point &p1,Point &p2) - 求直线和圆的交点，返回交点个数 gercircle() - 得到过a,b两点，半径为r1的两个圆 getcircle() - 得到与直线u相切，过点q,半径为r1的圆 getcircle() - 同时与直线u,v相切，半径为r1的圆 getcircle() - 同时与不相交圆cx,cy相切，半径为r1的圆 tangentline(Point q,Line &u,Line &v) - 过一点作圆的切线(先判断点和圆的关系) areacircle(circle v) - 求两圆相交的面积 areatriangle(Point a,Point b) - 求圆和三角形pab的相交面积 struct circle{ Point p;//圆心 double r;//半径 circle(){} circle(Point _p,double _r){ p = _p; r = _r; } circle(double x,double y,double _r){ p = Point(x,y); r = _r; } //`三角形的外接圆` //`需要Point的+ / rotate() 以及Line的crosspoint()` //`利用两条边的中垂线得到圆心` //`测试：UVA12304` circle(Point a,Point b,Point c){ Line u = Line((a+b)/2,((a+b)/2)+((b-a).rotleft())); Line v = Line((b+c)/2,((b+c)/2)+((c-b).rotleft())); p = u.crosspoint(v); r = p.distance(a); } //`三角形的内切圆` //`参数bool t没有作用，只是为了和上面外接圆函数区别` //`测试：UVA12304` circle(Point a,Point b,Point c,bool t){ Line u,v; double m = atan2(b.y-a.y,b.x-a.x), n = atan2(c.y-a.y,c.x-a.x); u.s = a; u.e = u.s + Point(cos((n+m)/2),sin((n+m)/2)); v.s = b; m = atan2(a.y-b.y,a.x-b.x) , n = atan2(c.y-b.y,c.x-b.x); v.e = v.s + Point(cos((n+m)/2),sin((n+m)/2)); p = u.crosspoint(v); r = Line(a,b).dispointtoseg(p); } //输入 void input(){ p.input(); scanf("%lf",&r); } //输出 void output(){ printf("%.2lf %.2lf %.2lf\n",p.x,p.y,r); } bool operator == (circle v){ return (p==v.p) && sgn(r-v.r)==0; } bool operator < (circle v)const{ return ((p<v.p)||((p==v.p)&&sgn(r-v.r)<0)); } //面积 double area(){ return pi*r*r; } //周长 double circumference(){ return 2*pi*r; } //`点和圆的关系` //`0 圆外` //`1 圆上` //`2 圆内` int relation(Point b){ double dst = b.distance(p); if(sgn(dst-r) < 0)return 2; else if(sgn(dst-r)==0)return 1; return 0; } //`线段和圆的关系` //`比较的是圆心到线段的距离和半径的关系` int relationseg(Line v){ double dst = v.dispointtoseg(p); if(sgn(dst-r) < 0)return 2; else if(sgn(dst-r) == 0)return 1; return 0; } //`直线和圆的关系` //`比较的是圆心到直线的距离和半径的关系` int relationline(Line v){ double dst = v.dispointtoline(p); if(sgn(dst-r) < 0)return 2; else if(sgn(dst-r) == 0)return 1; return 0; } //`两圆的关系` //`5 相离` //`4 外切` //`3 相交` //`2 内切` //`1 内含` //`需要Point的distance` //`测试：UVA12304` int relationcircle(circle v){ double d = p.distance(v.p); if(sgn(d-r-v.r) > 0) return 5; if(sgn(d-r-v.r) == 0) return 4; double l = fabs(r - v.r); if(sgn(d-r-v.r)<0 && sgn(d-l)>0) return 3; if(sgn(d-l)==0) return 2; if(sgn(d-l)<0) return 1; } //`求两个圆的交点，返回0表示没有交点，返回1是一个交点，2是两个交点` //`需要relationcircle` //`测试：UVA12304` int pointcrosscircle(circle v,Point &p1,Point &p2){ int rel = relationcircle(v); if(rel == 1 || rel == 5)return 0; double d = p.distance(v.p); double l = (d*d+r*r-v.r*v.r)/(2*d); double h = sqrt(r*r-l*l); Point tmp = p + (v.p-p).trunc(l); p1 = tmp + ((v.p-p).rotleft().trunc(h)); p2 = tmp + ((v.p-p).rotright().trunc(h)); if(rel == 2 || rel == 4) return 1; return 2; } //`求直线和圆的交点，返回交点个数` int pointcrossline(Line v,Point &p1,Point &p2){ if(!(*this).relationline(v))return 0; Point a = v.lineprog(p); double d = v.dispointtoline(p); d = sqrt(r*r-d*d); if(sgn(d) == 0){ p1 = a; p2 = a; return 1; } p1 = a + (v.e-v.s).trunc(d); p2 = a - (v.e-v.s).trunc(d); return 2; } //`得到过a,b两点，半径为r1的两个圆` int gercircle(Point a,Point b,double r1,circle &c1,circle &c2){ circle x(a,r1),y(b,r1); int t = x.pointcrosscircle(y,c1.p,c2.p); if(!t)return 0; c1.r = c2.r = r; return t; } //`得到与直线u相切，过点q,半径为r1的圆` //`测试：UVA12304` int getcircle(Line u,Point q,double r1,circle &c1,circle &c2){ double dis = u.dispointtoline(q); if(sgn(dis-r1*2)>0)return 0; if(sgn(dis) == 0){ c1.p = q + ((u.e-u.s).rotleft().trunc(r1)); c2.p = q + ((u.e-u.s).rotright().trunc(r1)); c1.r = c2.r = r1; return 2; } Line u1 = Line((u.s + (u.e-u.s).rotleft().trunc(r1)),(u.e + (u.e-u.s).rotleft().trunc(r1))); Line u2 = Line((u.s + (u.e-u.s).rotright().trunc(r1)),(u.e + (u.e-u.s).rotright().trunc(r1))); circle cc = circle(q,r1); Point p1,p2; if(!cc.pointcrossline(u1,p1,p2))cc.pointcrossline(u2,p1,p2); c1 = circle(p1,r1); if(p1 == p2){ c2 = c1; return 1; } c2 = circle(p2,r1); return 2; } //`同时与直线u,v相切，半径为r1的圆` //`测试：UVA12304` int getcircle(Line u,Line v,double r1,circle &c1,circle &c2,circle &c3,circle &c4){ if(u.parallel(v))return 0;//两直线平行 Line u1 = Line(u.s + (u.e-u.s).rotleft().trunc(r1),u.e + (u.e-u.s).rotleft().trunc(r1)); Line u2 = Line(u.s + (u.e-u.s).rotright().trunc(r1),u.e + (u.e-u.s).rotright().trunc(r1)); Line v1 = Line(v.s + (v.e-v.s).rotleft().trunc(r1),v.e + (v.e-v.s).rotleft().trunc(r1)); Line v2 = Line(v.s + (v.e-v.s).rotright().trunc(r1),v.e + (v.e-v.s).rotright().trunc(r1)); c1.r = c2.r = c3.r = c4.r = r1; c1.p = u1.crosspoint(v1); c2.p = u1.crosspoint(v2); c3.p = u2.crosspoint(v1); c4.p = u2.crosspoint(v2); return 4; } //`同时与不相交圆cx,cy相切，半径为r1的圆` //`测试：UVA12304` int getcircle(circle cx,circle cy,double r1,circle &c1,circle &c2){ circle x(cx.p,r1+cx.r),y(cy.p,r1+cy.r); int t = x.pointcrosscircle(y,c1.p,c2.p); if(!t)return 0; c1.r = c2.r = r1; return t; } //`过一点作圆的切线(先判断点和圆的关系)` //`测试：UVA12304` int tangentline(Point q,Line &u,Line &v){ int x = relation(q); if(x == 2)return 0; if(x == 1){ u = Line(q,q + (q-p).rotleft()); v = u; return 1; } double d = p.distance(q); double l = r*r/d; double h = sqrt(r*r-l*l); u = Line(q,p + ((q-p).trunc(l) + (q-p).rotleft().trunc(h))); v = Line(q,p + ((q-p).trunc(l) + (q-p).rotright().trunc(h))); return 2; } //`求两圆相交的面积` double areacircle(circle v){ int rel = relationcircle(v); if(rel >= 4)return 0.0; if(rel <= 2)return min(area(),v.area()); double d = p.distance(v.p); double hf = (r+v.r+d)/2.0; double ss = 2*sqrt(hf*(hf-r)*(hf-v.r)*(hf-d)); double a1 = acos((r*r+d*d-v.r*v.r)/(2.0*r*d)); a1 = a1*r*r; double a2 = acos((v.r*v.r+d*d-r*r)/(2.0*v.r*d)); a2 = a2*v.r*v.r; return a1+a2-ss; } //`求圆和三角形pab的相交面积` //`测试：POJ3675 HDU3982 HDU2892` double areatriangle(Point a,Point b){ if(sgn((p-a)^(p-b)) == 0)return 0.0; Point q[5]; int len = 0; q[len++] = a; Line l(a,b); Point p1,p2; if(pointcrossline(l,q[1],q[2])==2){ if(sgn((a-q[1])*(b-q[1]))<0)q[len++] = q[1]; if(sgn((a-q[2])*(b-q[2]))<0)q[len++] = q[2]; } q[len++] = b; if(len == 4 && sgn((q[0]-q[1])*(q[2]-q[1]))>0)swap(q[1],q[2]); double res = 0; for(int i = 0;i < len-1;i++){ if(relation(q[i])==0||relation(q[i+1])==0){ double arg = p.rad(q[i],q[i+1]); res += r*r*arg/2.0; } else{ res += fabs((q[i]-p)^(q[i+1]-p))/2.0; } } return res; } }; 第四部分—多边形 n,p line l for each side input(int _n) - 输入一个n边形 add(point q) - 在列表末尾添加一个点 getline() - populates line array cmp - comparision in convex_hull order norm() - 进行极角排序,首先需要找到最左下角的点,需要重载号好point的 < 操作符 getconvex(polygon &convex) - 得到凸包 注意如果有影响，要特判下所有点共点，或者共线的特殊情况 Graham(polygon &convex) - 得到凸包的另外一种方法 isconvex() - 判断是不是凸的 relationpoint(point q) - 判断点和任意多边形的关系 3 点上 2 边上 1 内部 0 外部 convexcut(line u,polygon &po) - 直线u切割凸多边形左侧，注意直线方向 gercircumference() - returns 周长 //测试 LightOJ1239 getarea() - returns 面积 getdir() - returns 方向：0表示顺时针, 1表示逆时针 getbarycentre() - returns 重心 areacircle(circle c) - 多边形和圆交的面积 relationcircle(circle c) - 多边形和圆关系 2 圆完全在多边形内 1 圆在多边形里面，碰到了多边形边界 0 其它 struct polygon{ int n; Point p[maxp]; Line l[maxp]; void input(int _n){ n = _n; for(int i = 0;i < n;i++) p[i].input(); } void add(Point q){ p[n++] = q; } void getline(){ for(int i = 0;i < n;i++){ l[i] = Line(p[i],p[(i+1)%n]); } } struct cmp{ Point p; cmp(const Point &p0){p = p0;} bool operator()(const Point &aa,const Point &bb){ Point a = aa, b = bb; int d = sgn((a-p)^(b-p)); if(d == 0){ return sgn(a.distance(p)-b.distance(p)) < 0; } return d > 0; } }; //`进行极角排序` //`首先需要找到最左下角的点` //`需要重载号好Point的 < 操作符(min函数要用) ` void norm(){ Point mi = p[0]; for(int i = 1;i < n;i++)mi = min(mi,p[i]); sort(p,p+n,cmp(mi)); } //`得到凸包` //`得到的凸包里面的点编号是0$\sim$n-1的` //`两种凸包的方法` //`注意如果有影响，要特判下所有点共点，或者共线的特殊情况` //`测试 LightOJ1203 LightOJ1239` void getconvex(polygon &convex){ sort(p,p+n); convex.n = n; for(int i = 0;i < min(n,2);i++){ convex.p[i] = p[i]; } if(convex.n == 2 && (convex.p[0] == convex.p[1]))convex.n--;//特判 if(n <= 2)return; int &top = convex.n; top = 1; for(int i = 2;i < n;i++){ while(top && sgn((convex.p[top]-p[i])^(convex.p[top-1]-p[i])) <= 0) top--; convex.p[++top] = p[i]; } int temp = top; convex.p[++top] = p[n-2]; for(int i = n-3;i >= 0;i--){ while(top != temp && sgn((convex.p[top]-p[i])^(convex.p[top-1]-p[i])) <= 0) top--; convex.p[++top] = p[i]; } if(convex.n == 2 && (convex.p[0] == convex.p[1]))convex.n--;//特判 convex.norm();//`原来得到的是顺时针的点，排序后逆时针` } //`得到凸包的另外一种方法` //`测试 LightOJ1203 LightOJ1239` void Graham(polygon &convex){ norm(); int &top = convex.n; top = 0; if(n == 1){ top = 1; convex.p[0] = p[0]; return; } if(n == 2){ top = 2; convex.p[0] = p[0]; convex.p[1] = p[1]; if(convex.p[0] == convex.p[1])top--; return; } convex.p[0] = p[0]; convex.p[1] = p[1]; top = 2; for(int i = 2;i < n;i++){ while( top > 1 && sgn((convex.p[top-1]-convex.p[top-2])^(p[i]-convex.p[top-2])) <= 0 ) top--; convex.p[top++] = p[i]; } if(convex.n == 2 && (convex.p[0] == convex.p[1]))convex.n--;//特判 } //`判断是不是凸的` bool isconvex() { bool s[3]; memset(s, false, sizeof(s)); for(int i = 0;i < n;i++){ int j = (i+1)%n; int k = (j+1)%n; s[sgn((p[j]-p[i])^(p[k]-p[i]))+1] = true; if(s[0] && s[2])return false; } return true; } //`判断点和任意多边形的关系` //` 3 点上` //` 2 边上` //` 1 内部` //` 0 外部` int relationpoint(Point q){ for(int i = 0;i < n;i++){ if(p[i] == q)return 3; } getline(); for(int i = 0;i < n;i++){ if(l[i].pointonseg(q))return 2; } int cnt = 0; for(int i = 0;i < n;i++){ int j = (i+1)%n; int k = sgn((q-p[j])^(p[i]-p[j])); int u = sgn(p[i].y-q.y); int v = sgn(p[j].y-q.y); if(k > 0 && u < 0 && v >= 0)cnt++; if(k < 0 && v < 0 && u >= 0)cnt--; } return cnt != 0; } //`直线u切割凸多边形左侧` //`注意直线方向` //`测试：HDU3982` void convexcut(Line u,polygon &po){ int &top = po.n;//注意引用 top = 0; for(int i = 0;i < n;i++){ int d1 = sgn((u.e-u.s)^(p[i]-u.s)); int d2 = sgn((u.e-u.s)^(p[(i+1)%n]-u.s)); if(d1 >= 0)po.p[top++] = p[i]; if(d1*d2 < 0)po.p[top++] = u.crosspoint(Line(p[i],p[(i+1)%n])); } } //`得到周长` //`测试 LightOJ1239` double getcircumference(){ double sum = 0; for(int i = 0;i < n;i++){ sum += p[i].distance(p[(i+1)%n]); } return sum; } //`得到面积` double getarea(){ double sum = 0; for(int i = 0;i < n;i++){ sum += (p[i]^p[(i+1)%n]); } return fabs(sum)/2; } //`得到方向` //` 1 表示逆时针，0表示顺时针` bool getdir(){ double sum = 0; for(int i = 0;i < n;i++) sum += (p[i]^p[(i+1)%n]); if(sgn(sum) > 0)return 1; return 0; } //`得到重心` Point getbarycentre(){ Point ret(0,0); double area = 0; for(int i = 1;i < n-1;i++){ double tmp = (p[i]-p[0])^(p[i+1]-p[0]); if(sgn(tmp) == 0)continue; area += tmp; ret.x += (p[0].x+p[i].x+p[i+1].x)/3*tmp; ret.y += (p[0].y+p[i].y+p[i+1].y)/3*tmp; } if(sgn(area)) ret = ret/area; return ret; } //`多边形和圆交的面积` //`测试：POJ3675 HDU3982 HDU2892` double areacircle(circle c){ double ans = 0; for(int i = 0;i < n;i++){ int j = (i+1)%n; if(sgn( (p[j]-c.p)^(p[i]-c.p) ) >= 0) ans += c.areatriangle(p[i],p[j]); else ans -= c.areatriangle(p[i],p[j]); } return fabs(ans); } //`多边形和圆关系` //` 2 圆完全在多边形内` //` 1 圆在多边形里面，碰到了多边形边界` //` 0 其它` int relationcircle(circle c){ getline(); int x = 2; if(relationpoint(c.p) != 1)return 0;//圆心不在内部 for(int i = 0;i < n;i++){ if(c.relationseg(l[i])==2)return 0; if(c.relationseg(l[i])==1)x = 1; } return x; } }; 其他 //`AB X AC` double cross(Point A,Point B,Point C){ return (B-A)^(C-A); } //`AB*AC` double dot(Point A,Point B,Point C){ return (B-A)*(C-A); } 最小矩形面积覆盖 //` A 必须是凸包(而且是逆时针顺序)` //` 测试 UVA 10173` double minRectangleCover(polygon A){ //`要特判A.n < 3的情况` if(A.n < 3)return 0.0; A.p[A.n] = A.p[0]; double ans = -1; int r = 1, p = 1, q; for(int i = 0;i < A.n;i++){ //`卡出离边A.p[i] - A.p[i+1]最远的点` while( sgn( cross(A.p[i],A.p[i+1],A.p[r+1]) - cross(A.p[i],A.p[i+1],A.p[r]) ) >= 0 ) r = (r+1)%A.n; //`卡出A.p[i] - A.p[i+1]方向上正向n最远的点` while(sgn( dot(A.p[i],A.p[i+1],A.p[p+1]) - dot(A.p[i],A.p[i+1],A.p[p]) ) >= 0 ) p = (p+1)%A.n; if(i == 0)q = p; //`卡出A.p[i] - A.p[i+1]方向上负向最远的点` while(sgn(dot(A.p[i],A.p[i+1],A.p[q+1]) - dot(A.p[i],A.p[i+1],A.p[q])) <= 0) q = (q+1)%A.n; double d = (A.p[i] - A.p[i+1]).len2(); double tmp = cross(A.p[i],A.p[i+1],A.p[r]) * (dot(A.p[i],A.p[i+1],A.p[p]) - dot(A.p[i],A.p[i+1],A.p[q]))/d; if(ans < 0 || ans > tmp)ans = tmp; } return ans; } 直线切凸多边形 //`多边形是逆时针的，在q1q2的左侧` //`测试:HDU3982` vector<Point> convexCut(const vector<Point> &ps,Point q1,Point q2){ vector<Point>qs; int n = ps.size(); for(int i = 0;i < n;i++){ Point p1 = ps[i], p2 = ps[(i+1)%n]; int d1 = sgn((q2-q1)^(p1-q1)), d2 = sgn((q2-q1)^(p2-q1)); if(d1 >= 0) qs.push_back(p1); if(d1 * d2 < 0) qs.push_back(Line(p1,p2).crosspoint(Line(q1,q2))); } return qs; } 半平面交 //`测试 POJ3335 POJ1474 POJ1279` //*************************** struct halfplane:public Line{ double angle; halfplane(){} //`表示向量s->e逆时针(左侧)的半平面` halfplane(Point _s,Point _e){ s = _s; e = _e; } halfplane(Line v){ s = v.s; e = v.e; } void calcangle(){ angle = atan2(e.y-s.y,e.x-s.x); } bool operator <(const halfplane &b)const{ return angle < b.angle; } }; struct halfplanes{ int n; halfplane hp[maxp]; Point p[maxp]; int que[maxp]; int st,ed; void push(halfplane tmp){ hp[n++] = tmp; } //去重 void unique(){ int m = 1; for(int i = 1;i < n;i++){ if(sgn(hp[i].angle-hp[i-1].angle) != 0) hp[m++] = hp[i]; else if(sgn( (hp[m-1].e-hp[m-1].s)^(hp[i].s-hp[m-1].s) ) > 0) hp[m-1] = hp[i]; } n = m; } bool halfplaneinsert(){ for(int i = 0;i < n;i++)hp[i].calcangle(); sort(hp,hp+n); unique(); que[st=0] = 0; que[ed=1] = 1; p[1] = hp[0].crosspoint(hp[1]); for(int i = 2;i < n;i++){ while(st<ed && sgn((hp[i].e-hp[i].s)^(p[ed]-hp[i].s))<0)ed--; while(st<ed && sgn((hp[i].e-hp[i].s)^(p[st+1]-hp[i].s))<0)st++; que[++ed] = i; if(hp[i].parallel(hp[que[ed-1]]))return false; p[ed]=hp[i].crosspoint(hp[que[ed-1]]); } while(st<ed && sgn((hp[que[st]].e-hp[que[st]].s)^(p[ed]-hp[que[st]].s))<0)ed--; while(st<ed && sgn((hp[que[ed]].e-hp[que[ed]].s)^(p[st+1]-hp[que[ed]].s))<0)st++; if(st+1>=ed)return false; return true; } //`得到最后半平面交得到的凸多边形` //`需要先调用halfplaneinsert() 且返回true` void getconvex(polygon &con){ p[st] = hp[que[st]].crosspoint(hp[que[ed]]); con.n = ed-st+1; for(int j = st,i = 0;j <= ed;i++,j++) con.p[i] = p[j]; } }; const int maxn = 1010; struct circles{ circle c[maxn]; double ans[maxn];//`ans[i]表示被覆盖了i次的面积` double pre[maxn]; int n; circles(){} void add(circle cc){ c[n++] = cc; } //`x包含在y中` bool inner(circle x,circle y){ if(x.relationcircle(y) != 1)return 0; return sgn(x.r-y.r)<=0?1:0; } //圆的面积并去掉内含的圆 void init_or(){ bool mark[maxn] = {0}; int i,j,k=0; for(i = 0;i < n;i++){ for(j = 0;j < n;j++) if(i != j && !mark[j]){ if( (c[i]==c[j])||inner(c[i],c[j]) )break; } if(j < n)mark[i] = 1; } for(i = 0;i < n;i++) if(!mark[i]) c[k++] = c[i]; n = k; } //`圆的面积交去掉内含的圆` void init_add(){ int i,j,k; bool mark[maxn] = {0}; for(i = 0;i < n;i++){ for(j = 0;j < n;j++) if(i != j && !mark[j]){ if( (c[i]==c[j])||inner(c[j],c[i]) )break; } if(j < n)mark[i] = 1; } for(i = 0;i < n;i++) if(!mark[i]) c[k++] = c[i]; n = k; } //`半径为r的圆，弧度为th对应的弓形的面积` double areaarc(double th,double r){ return 0.5*r*r*(th-sin(th)); } //`测试SPOJVCIRCLES SPOJCIRUT` //`SPOJVCIRCLES求n个圆并的面积，需要加上init\_or()去掉重复圆（否则WA）` //`SPOJCIRUT 是求被覆盖k次的面积，不能加init\_or()` //`对于求覆盖多少次面积的问题，不能解决相同圆，而且不能init\_or()` //`求多圆面积并，需要init\_or,其中一个目的就是去掉相同圆` void getarea(){ memset(ans,0,sizeof(ans)); vector<pair<double,int> >v; for(int i = 0;i < n;i++){ v.clear(); v.push_back(make_pair(-pi,1)); v.push_back(make_pair(pi,-1)); for(int j = 0;j < n;j++) if(i != j){ Point q = (c[j].p - c[i].p); double ab = q.len(),ac = c[i].r, bc = c[j].r; if(sgn(ab+ac-bc)<=0){ v.push_back(make_pair(-pi,1)); v.push_back(make_pair(pi,-1)); continue; } if(sgn(ab+bc-ac)<=0)continue; if(sgn(ab-ac-bc)>0)continue; double th = atan2(q.y,q.x), fai = acos((ac*ac+ab*ab-bc*bc)/(2.0*ac*ab)); double a0 = th-fai; if(sgn(a0+pi)<0)a0+=2*pi; double a1 = th+fai; if(sgn(a1-pi)>0)a1-=2*pi; if(sgn(a0-a1)>0){ v.push_back(make_pair(a0,1)); v.push_back(make_pair(pi,-1)); v.push_back(make_pair(-pi,1)); v.push_back(make_pair(a1,-1)); } else{ v.push_back(make_pair(a0,1)); v.push_back(make_pair(a1,-1)); } } sort(v.begin(),v.end()); int cur = 0; for(int j = 0;j < v.size();j++){ if(cur && sgn(v[j].first-pre[cur])){ ans[cur] += areaarc(v[j].first-pre[cur],c[i].r); ans[cur] += 0.5*(Point(c[i].p.x+c[i].r*cos(pre[cur]),c[i].p.y+c[i].r*sin(pre[cur]))^Point(c[i].p.x+c[i].r*cos(v[j].first),c[i].p.y+c[i].r*sin(v[j].first))); } cur += v[j].second; pre[cur] = v[j].first; } } for(int i = 1;i < n;i++) ans[i] -= ans[i+1]; } }; 三维几何 const double eps = 1e-8; int sgn(double x){ if(fabs(x) < eps)return 0; if(x < 0)return -1; else return 1; } struct Point3{ double x,y,z; Point3(double _x = 0,double _y = 0,double _z = 0){ x = _x; y = _y; z = _z; } void input(){ scanf("%lf%lf%lf",&x,&y,&z); } void output(){ scanf("%.2lf %.2lf %.2lf\n",x,y,z); } bool operator ==(const Point3 &b)const{ return sgn(x-b.x) == 0 && sgn(y-b.y) == 0 && sgn(z-b.z) == 0; } bool operator <(const Point3 &b)const{ return sgn(x-b.x)==0?(sgn(y-b.y)==0?sgn(z-b.z)<0:y<b.y):x<b.x; } double len(){ return sqrt(x*x+y*y+z*z); } double len2(){ return x*x+y*y+z*z; } double distance(const Point3 &b)const{ return sqrt((x-b.x)*(x-b.x)+(y-b.y)*(y-b.y)+(z-b.z)*(z-b.z)); } Point3 operator -(const Point3 &b)const{ return Point3(x-b.x,y-b.y,z-b.z); } Point3 operator +(const Point3 &b)const{ return Point3(x+b.x,y+b.y,z+b.z); } Point3 operator *(const double &k)const{ return Point3(x*k,y*k,z*k); } Point3 operator /(const double &k)const{ return Point3(x/k,y/k,z/k); } //点乘 double operator *(const Point3 &b)const{ return x*b.x+y*b.y+z*b.z; } //叉乘 Point3 operator ^(const Point3 &b)const{ return Point3(y*b.z-z*b.y,z*b.x-x*b.z,x*b.y-y*b.x); } double rad(Point3 a,Point3 b){ Point3 p = (*this); return acos( ( (a-p)*(b-p) )/ (a.distance(p)*b.distance(p)) ); } //变换长度 Point3 trunc(double r){ double l = len(); if(!sgn(l))return *this; r /= l; return Point3(x*r,y*r,z*r); } }; struct Line3 { Point3 s,e; Line3(){} Line3(Point3 _s,Point3 _e) { s = _s; e = _e; } bool operator ==(const Line3 v) { return (s==v.s)&&(e==v.e); } void input() { s.input(); e.input(); } double length() { return s.distance(e); } //点到直线距离 double dispointtoline(Point3 p) { return ((e-s)^(p-s)).len()/s.distance(e); } //点到线段距离 double dispointtoseg(Point3 p) { if(sgn((p-s)*(e-s)) < 0 || sgn((p-e)*(s-e)) < 0) return min(p.distance(s),e.distance(p)); return dispointtoline(p); } //`返回点p在直线上的投影` Point3 lineprog(Point3 p) { return s + ( ((e-s)*((e-s)*(p-s)))/((e-s).len2()) ); } //`p绕此向量逆时针arg角度` Point3 rotate(Point3 p,double ang) { if(sgn(((s-p)^(e-p)).len()) == 0)return p; Point3 f1 = (e-s)^(p-s); Point3 f2 = (e-s)^(f1); double len = ((s-p)^(e-p)).len()/s.distance(e); f1 = f1.trunc(len); f2 = f2.trunc(len); Point3 h = p+f2; Point3 pp = h+f1; return h + ((p-h)*cos(ang)) + ((pp-h)*sin(ang)); } //`点在直线上` bool pointonseg(Point3 p) { return sgn( ((s-p)^(e-p)).len() ) == 0 && sgn((s-p)*(e-p)) == 0; } }; struct Plane { Point3 a,b,c,o;//`平面上的三个点，以及法向量` Plane(){} Plane(Point3 _a,Point3 _b,Point3 _c) { a = _a; b = _b; c = _c; o = pvec(); } Point3 pvec() { return (b-a)^(c-a); } //`ax+by+cz+d = 0` Plane(double _a,double _b,double _c,double _d) { o = Point3(_a,_b,_c); if(sgn(_a) != 0) a = Point3((-_d-_c-_b)/_a,1,1); else if(sgn(_b) != 0) a = Point3(1,(-_d-_c-_a)/_b,1); else if(sgn(_c) != 0) a = Point3(1,1,(-_d-_a-_b)/_c); } //`点在平面上的判断` bool pointonplane(Point3 p) { return sgn((p-a)*o) == 0; } //`两平面夹角` double angleplane(Plane f) { return acos(o*f.o)/(o.len()*f.o.len()); } //`平面和直线的交点，返回值是交点个数` int crossline(Line3 u,Point3 &p) { double x = o*(u.e-a); double y = o*(u.s-a); double d = x-y; if(sgn(d) == 0)return 0; p = ((u.s*x)-(u.e*y))/d; return 1; } //`点到平面最近点(也就是投影)` Point3 pointtoplane(Point3 p) { Line3 u = Line3(p,p+o); crossline(u,p); return p; } //`平面和平面的交线` int crossplane(Plane f,Line3 &u) { Point3 oo = o^f.o; Point3 v = o^oo; double d = fabs(f.o*v); if(sgn(d) == 0)return 0; Point3 q = a + (v*(f.o*(f.a-a))/d); u = Line3(q,q+oo); return 1; } }; 平面最近点对 const int MAXN = 100010; const double eps = 1e-8; const double INF = 1e20; struct Point{ double x,y; void input(){ scanf("%lf%lf",&x,&y); } }; double dist(Point a,Point b){ return sqrt((a.x-b.x)*(a.x-b.x) + (a.y-b.y)*(a.y-b.y)); } Point p[MAXN]; Point tmpt[MAXN]; bool cmpx(Point a,Point b){ return a.x < b.x || (a.x == b.x && a.y < b.y); } bool cmpy(Point a,Point b){ return a.y < b.y || (a.y == b.y && a.x < b.x); } double Closest_Pair(int left,int right){ double d = INF; if(left == right)return d; if(left+1 == right)return dist(p[left],p[right]); int mid = (left+right)/2; double d1 = Closest_Pair(left,mid); double d2 = Closest_Pair(mid+1,right); d = min(d1,d2); int cnt = 0; for(int i = left;i <= right;i++){ if(fabs(p[mid].x - p[i].x) <= d) tmpt[cnt++] = p[i]; } sort(tmpt,tmpt+cnt,cmpy); for(int i = 0;i < cnt;i++){ for(int j = i+1;j < cnt && tmpt[j].y - tmpt[i].y < d;j++) d = min(d,dist(tmpt[i],tmpt[j])); } return d; } int main(){ int n; while(scanf("%d",&n) == 1 && n){ for(int i = 0;i < n;i++)p[i].input(); sort(p,p+n,cmpx); printf("%.2lf\n",Closest_Pair(0,n-1)); } return 0; } 三维凸包   例题：HDU - 4273 const double eps = 1e-8; const int MAXN = 550; int sgn(double x){ if(fabs(x) < eps)return 0; if(x < 0)return -1; else return 1; } struct Point3{ double x,y,z; Point3(double _x = 0, double _y = 0, double _z = 0){ x = _x; y = _y; z = _z; } void input(){ scanf("%lf%lf%lf",&x,&y,&z); } bool operator ==(const Point3 &b)const{ return sgn(x-b.x) == 0 && sgn(y-b.y) == 0 && sgn(z-b.z) == 0; } double len(){ return sqrt(x*x+y*y+z*z); } double len2(){ return x*x+y*y+z*z; } double distance(const Point3 &b)const{ return sqrt((x-b.x)*(x-b.x)+(y-b.y)*(y-b.y)+(z-b.z)*(z-b.z)); } Point3 operator -(const Point3 &b)const{ return Point3(x-b.x,y-b.y,z-b.z); } Point3 operator +(const Point3 &b)const{ return Point3(x+b.x,y+b.y,z+b.z); } Point3 operator *(const double &k)const{ return Point3(x*k,y*k,z*k); } Point3 operator /(const double &k)const{ return Point3(x/k,y/k,z/k); } //点乘 double operator *(const Point3 &b)const{ return x*b.x + y*b.y + z*b.z; } //叉乘 Point3 operator ^(const Point3 &b)const{ return Point3(y*b.z-z*b.y,z*b.x-x*b.z,x*b.y-y*b.x); } }; struct CH3D{ struct face{ //表示凸包一个面上的三个点的编号 int a,b,c; //表示该面是否属于最终的凸包上的面 bool ok; }; //初始顶点数 int n; Point3 P[MAXN]; //凸包表面的三角形数 int num; //凸包表面的三角形 face F[8*MAXN]; int g[MAXN][MAXN]; //叉乘 Point3 cross(const Point3 &a,const Point3 &b,const Point3 &c){ return (b-a)^(c-a); } //`三角形面积*2` double area(Point3 a,Point3 b,Point3 c){ return ((b-a)^(c-a)).len(); } //`四面体有向面积*6` double volume(Point3 a,Point3 b,Point3 c,Point3 d){ return ((b-a)^(c-a))*(d-a); } //`正：点在面同向` double dblcmp(Point3 &p,face &f){ Point3 p1 = P[f.b] - P[f.a]; Point3 p2 = P[f.c] - P[f.a]; Point3 p3 = p - P[f.a]; return (p1^p2)*p3; } void deal(int p,int a,int b){ int f = g[a][b]; face add; if(F[f].ok){ if(dblcmp(P[p],F[f]) > eps) dfs(p,f); else { add.a = b; add.b = a; add.c = p; add.ok = true; g[p][b] = g[a][p] = g[b][a] = num; F[num++] = add; } } } //递归搜索所有应该从凸包内删除的面 void dfs(int p,int now){ F[now].ok = false; deal(p,F[now].b,F[now].a); deal(p,F[now].c,F[now].b); deal(p,F[now].a,F[now].c); } bool same(int s,int t){ Point3 &a = P[F[s].a]; Point3 &b = P[F[s].b]; Point3 &c = P[F[s].c]; return fabs(volume(a,b,c,P[F[t].a])) < eps && fabs(volume(a,b,c,P[F[t].b])) < eps && fabs(volume(a,b,c,P[F[t].c])) < eps; } //构建三维凸包 void create(){ num = 0; face add; //*********************************** //此段是为了保证前四个点不共面 bool flag = true; for(int i = 1;i < n;i++){ if(!(P[0] == P[i])){ swap(P[1],P[i]); flag = false; break; } } if(flag)return; flag = true; for(int i = 2;i < n;i++){ if( ((P[1]-P[0])^(P[i]-P[0])).len() > eps ){ swap(P[2],P[i]); flag = false; break; } } if(flag)return; flag = true; for(int i = 3;i < n;i++){ if(fabs( ((P[1]-P[0])^(P[2]-P[0]))*(P[i]-P[0]) ) > eps){ swap(P[3],P[i]); flag = false; break; } } if(flag)return; //********************************** for(int i = 0;i < 4;i++){ add.a = (i+1)%4; add.b = (i+2)%4; add.c = (i+3)%4; add.ok = true; if(dblcmp(P[i],add) > 0)swap(add.b,add.c); g[add.a][add.b] = g[add.b][add.c] = g[add.c][add.a] = num; F[num++] = add; } for(int i = 4;i < n;i++) for(int j = 0;j < num;j++) if(F[j].ok && dblcmp(P[i],F[j]) > eps){ dfs(i,j); break; } int tmp = num; num = 0; for(int i = 0;i < tmp;i++) if(F[i].ok) F[num++] = F[i]; } //表面积 //`测试：HDU3528` double area(){ double res = 0; if(n == 3){ Point3 p = cross(P[0],P[1],P[2]); return p.len()/2; } for(int i = 0;i < num;i++) res += area(P[F[i].a],P[F[i].b],P[F[i].c]); return res/2.0; } double volume(){ double res = 0; Point3 tmp = Point3(0,0,0); for(int i = 0;i < num;i++) res += volume(tmp,P[F[i].a],P[F[i].b],P[F[i].c]); return fabs(res/6); } //表面三角形个数 int triangle(){ return num; } //表面多边形个数 //`测试：HDU3662` int polygon(){ int res = 0; for(int i = 0;i < num;i++){ bool flag = true; for(int j = 0;j < i;j++) if(same(i,j)){ flag = 0; break; } res += flag; } return res; } //重心 //`测试：HDU4273` Point3 barycenter(){ Point3 ans = Point3(0,0,0); Point3 o = Point3(0,0,0); double all = 0; for(int i = 0;i < num;i++){ double vol = volume(o,P[F[i].a],P[F[i].b],P[F[i].c]); ans = ans + (((o+P[F[i].a]+P[F[i].b]+P[F[i].c])/4.0)*vol); all += vol; } ans = ans/all; return ans; } //点到面的距离 //`测试：HDU4273` double ptoface(Point3 p,int i){ double tmp1 = fabs(volume(P[F[i].a],P[F[i].b],P[F[i].c],p)); double tmp2 = ((P[F[i].b]-P[F[i].a])^(P[F[i].c]-P[F[i].a])).len(); return tmp1/tmp2; } }; CH3D hull; int main() { while(scanf("%d",&hull.n) == 1){ for(int i = 0;i < hull.n;i++)hull.P[i].input(); hull.create(); Point3 p = hull.barycenter(); double ans = 1e20; for(int i = 0;i < hull.num;i++) ans = min(ans,hull.ptoface(p,i)); printf("%.3lf\n",ans); } return 0; }

字符串算法总结整理 以前学习的做的笔记，现在整理一下发出来 字符串哈希 字典树 KMP next数组 板子 扩展KMP Manacher 字符串哈希 常用BKDRHash unsigned int BKDRHash(char *str) { unsigned int seed = 31, key = 0; while (*str) key = key * seed + (*str++); return key & 0x7fffffff; } ull base[maxn]; ull Hash[maxn]; base[0] = 1; for (int i = 1; i <= len; i++) base[i] = base[i - 1] * seed; Hash[0] = s[0] - 'a' + 1; for (int i = 1; i < len; i++) Hash[i] = Hash[i - 1] * seed + s[i] - 'a' + 1; 字典树 map使用:略 压缩空间字典树 int trie[maxn][26]; int num[maxn] = {0}; int pos = 1; void insert(char str[]) { int p = 0; for (int i = 0; str[i]; i++) { int n = str[i] - 'a'; if (trie[p][n] == 0) trie[p][n] = pos++; p = trie[p][n]; num[p]++; } } int find(char str[]) { int p = 0; for (int i = 0; str[i]; i++) { int n = str[i] - 'a'; if (trie[p][n] == 0) return 0; p = trie[p][n]; } } KMP KMP解决的是字符串匹配问题： 一个文本串S，一个模式串P，求P在S中出现的次数 next数组 next[i]为满足P[i-z ~ i-1] = P[0 ~ z-1] 的最大z值 next数组的含义：当前字符以前的字符串中，有多大长度的相同前缀和后缀 P[i-next[i]...i-1] = P[0 ... next[i] - 1] P[] A B C D A B D   next[] -1 0 0 0 0 1 2 0 void kmp_pre(char s[], int m, int next[]) { int i = 0, j = next[0] = -1; while (i < m) { while (j != -1 && x[i] != x[j]) j = next[j]; next[++i] = ++j; } } next数组的深入理解和应用： 沿着next数组一直往前后，得到的是所有前缀，也是当前后缀的字符串。 ABCABCABCAB -> ABCABCAB -> ABCAB -> AB P[] A B C A B C A B C A B   next[] -1 0 0 0 1 2 3 4 5 6 7 8 周期性字符串：n % (n-next[n]) == 0，循环节长度为n-nxet[n] next[]往前跳的步长是一样的，除了最后一次 板子 x[]为模式串, y[]为主串 int kmp_count(char x[], int m, char y[], int n) { int i = 0, j = 0, ans = 0; kmp_pre(x, m, next); while (i < n) { while (j != -1 && y[i] != x[j]) j = next[j]; i++; j++; if (j >= m) { ans++; j = next[j]; } } return ans; } 扩展KMP next[i]: x[i..m-1]与x[0...m-1]的最长公共前缀 extend[i]: y[i..n-1]与y[0///m-1]的最长公共前缀 例: y = "aaaaabbb" x = "aaaaac" next[] = {6, 4, 3, 2, 1, 0} extend[] = {5, 4, 3, 2, 1, 0, 0, 0} Manacher 求解问题：给定一个字符串s，找到s中最长的回文子字符串。 引入： 原串：abaaba 修改后的子串：$#a#b#a#a#b#a# 统一奇数长度和偶数长度的回文串，给所有字串一个中心 然后计算回文半径 i 0 1 2 3 4 5 6 7 8 9 10 11 12 13 修改后Ma[] $ # a # b # a # a # b # a # 回文半径Mp[] 1 1 2 1 4 1 2 7 2 1 4 1 2 1 Ma[i-Mp[i] ~ i+Mp[i]]是回文串 $符号可简化边界判断 算法步骤: 预处理得到字符串Ma 充分利用前面的信息得到数组Mp 时间复杂度和空间复杂度都是O(n) void Manacher(char s[], int len) { int l = 0; Ma[l++] = '$'; Ma[l++] = '#'; for (int i = 0; i < l; i++) { Ma[l++] = s[i]; Ma[l++] = '#'; } Ma[l] = 0; int mx = 0, id = 0; for (int i = 0; i < l; i++) { Mp[i] = mx > i ? min(Mp[2*id-i], mx - i) : 1; while (Ma[i+Mp[i]] == Ma[i-Mp[i]]) Mp[i]++; if (i = Mp[i] > mx) { mx = i + Mp[i]; id = i; } } }

数论基础算法整理总结 以前学习的做的笔记，现在整理一下发出来 一些概念 算术基本定理 费马定理 其他 快速幂 欧几里得求最大公约数 扩展欧几里得 中国剩余定理 欧拉函数和欧拉定理 线性基 定义 性质 插入 合并 查询 k小值 BM杜教筛 一些概念 算术基本定理 一个大于1的正整数 $N$，如果它的标准分解式为：$N = P^1_{a1} P^2_{a2}…P^n_{an}$，那么它的正因数个数为 $ \sigma_{0}(N) = (1+a_{1})(1+a_{2})...(1+a_{n}) $. 它的全体正因数之和为 $ \sigma_{1}(N) = (1+p_{1}+ p_{1}^{2}+...+ p_{1}^{a1})(1+p_{2}+ p_{2}^{2}+...+ p_{2}^{a2})...(1+p_{n}+ p_{n}^{2}+...+ p_{n}^{an}) $ 当 $ \sigma_{0}(N) = 2n $ 时就称 $N$ 为完全数。 是否存在奇完全数，是一个至今未解决之猜想。 利用算术基本定理可以重新定义整数 $a$ 和 $b$ 的最大公因子 $(a, b)$ 和最小公倍数 $[a, b]$，并证明 $ab = (a, b) × [a, b]$。 此外还可证明根号2是无理数等等。 证明素数个数无限。 费马定理 费马大定理：当 $n ≥ 2$ 时，关于$ x, y, z $的方程 $x^n + y^n = z^n$ 没有正整数解。 费马小定理： 假如 $p$ 是质数，且 $gcd(a, p)=1$，那么 $a^{p-1} ≡ 1（mod p）$，即：假如 $a$ 是整数，$p$ 是质数，且 $a, p$ 互质(即两者只有一个公约数1)，那么 $a$ 的 $(p-1)$ 次方除以 $p$ 的余数恒等于1。 其他 调和级数(即f(n))至今没有一个完全正确的公式，但欧拉给出过一个近似公式：(n很大时，n小的时候暴力) $f(n)≈ln(n)+C+1/2*n ( f(n)= 1/1 + 1/2 + 1/3 + … + 1/n )$ 欧拉常数值：$C≈0.57721566490153286060651209$ c++ math库中，log即为ln。 快速幂 ll pow(ll a, ll b) { ll result = 1; ll base = a; while (b) { if (b & 1) result *= base; base *= base; b >>= 1; } return result; } ll powmod(ll a, ll b, ll c) { //取模 ll result = 1; ll base = a % c; while (b) { if (b & 1) result = (result * base) % c; base = (base * base) % c; b >>= 1; } return result; } 欧几里得求最大公约数 ll gcd(ll a, ll b) { if (!b) return a; return gcd(b, a % b); } 扩展欧几里得 long long extend_gcd(long long a, long long b, long long &x, long long &y) { //求 ax + by = gcd(a, b) 的解 if (a == 0 && b == 0) return -1; //最大公约数不存在 if (b == 0) { x = 1, y = 0; return a; } long long r = extend_gcd(b, a % b, y, x); y -= x * (a / b); return r; } long long mod_inv(long long a, long long n) { //求逆元 long long x, y; long long d = extend_gcd(a, n, x, y); if (d == 1) return (x % n + n) % n; else return -1; } 中国剩余定理 //n个方程：x=a[i](mod m[i]) (0<=i<n) ll china(int n, ll *a, ll *m){ ll M = 1, ret = 0; for(int i = 0; i < n; i ++) M *= m[i]; for(int i = 0; i < n; i ++){ ll w = M / m[i]; ret = (ret + w * mod_inv(w, m[i]) * a[i]) % M; } return (ret + M) % M; } 欧拉函数和欧拉定理 定义：对正整数 n, 欧拉函数是小于 n 的正整数中与 n 互质的数的数目（φ(1)=1） 性质 : 对于质数 p, φ(p)=p−1 若p为质数, n=pk, 则 φ(n)=pk-pk-1 欧拉函数是积性函数, 但不是完全积性函数。若 m, n 互质, 则 φ(m∗n)=φ(m)∗φ(n) 特殊的, 当 m=2,n 为奇数时，φ(2*n)=φ(n) 当 n>2 时, φ(n)是偶数 小于 n 的数中, 与 n 互质的数的总和为：φ(n) * n / 2 (n>1) n=∑d∣nφ(d)，即n的因数 (包括 1 和它自己) 的欧拉函数之和等于 n (1) 求单个变量的欧拉函数 long long eular(long long x) { int res = x, a = x; for (int i = 2; i * i <= a; i++) { if (a % i == 0) { res = res / i * (i - 1); //res -= res/i; while (a % i == 0) a /= i; } } if(a > 1) res = res / a * (a - 1); //res -= res/a; return res; } (2) 欧拉表 int euler[Max]; void euler(void) { euler[1] = 1; for (int i = 2; i < Max; i++) euler[i] = i; for (int i = 2; i < Max; i++) if (euler[i] == i) for (int j = i; j < Max; j += i) euler[j] = euler[j] / i * (i - 1); //先进行除法是为了防止中间数据的溢出 } (3) 线性筛 (同时得到欧拉函数和素数表) const int MAXN = 10000000; bool check[MAXN + 10]; int phi[MAXN + 10]; int prime[MAXN + 10]; int tot; //素数的个数 void phi_and_prime_table(int N) { memset(check, false, sizeof(check)); phi[1] = 1; tot = 0; for (int i = 2; i <= N; i++) { if (!check[i]) { prime[tot++] = i; phi[i] = i - 1; } for (int j = 0; j < tot; j++) { if (i * prime[j] > N) break; check[i * prime[j]] = true; if (i % prime[j] == 0) { phi[i * prime[j]] = phi[i] * prime[j]; break; } else { phi[i * prime[j]] = phi[i] * (prime[j] - 1); } } } } 线性基 定义 设数集T的值域范围为[1, 2n−1]。 T的线性基是T的一个子集A={a1, a2, a3, ..., an}。 A中元素互相xor所形成的异或集合，等价于原数集T的元素互相xor形成的异或集合。 可以理解为将原数集进行了压缩。 性质 设线性基的异或集合中不存在0。 线性基的异或集合中每个元素的异或方案唯一，其实这个跟性质1是等价的。 线性基二进制最高位互不相同。 如果线性基是满的，它的异或集合为[1, 2n−1]。 线性基中元素互相异或，异或集合不变。 插入 如果向线性基中插入数x，从高位到低位扫描它为1的二进制位。 扫描到第i时，如果ai不存在，就令ai = x，否则x = x⊗ai。 x的结局是，要么被扔进线性基，要么经过一系列操作过后，变成了0。 bool insert(long long val) { for (int i = 60; i >= 0; i--) if (val & (1LL << i)) { if (!a[i]) { a[i] = val; break; } val ^= a[i]; } return val > 0; } 合并 将一个线性基暴力插入另一个线性基即可 L_B merge(const L_B& n1, const L_B& n2) { L_B ret = n1; for (int i = 0; i <= 60; i++) if (n2.d[i]) ret.insert(n2.d[i]); return ret; } 查询   如果要查询x是否存于异或集合中。   从高位到低位扫描x的为1的二进制位。   扫描到第i位的时候x=x⊗ai   如果中途x变为了0，那么表示x存于线性基的异或集合中。   从高位到低位扫描线性基。   如果异或后可以使得答案变大，就异或到答案中去。 long long query_max() { long long ret = 0; for (int i = 60; i >= 0; i--) if ((ret ^ d[i]) > ret) ret ^= d[i]; return ret; } 最大值   最小值即为最低位上的线性基 long long query_min() { for (int i = 0; i <= 60; i++) if (d[i]) return d[i]; return 0; } k小值 根据性质3。 我们要将线性基改造成每一位相互独立。 具体操作就是如果i < j，aj的第i位是1，就将aj异或上ai。 经过一系列操作之后，对于二进制的某一位i。只有ai的这一位是1，其他都是0。 所以查询的时候将k二进制拆分，对于1的位，就异或上对应的线性基。 最终得出的答案就是k小值。 void rebuild() { for (int i = 60; i >= 0; i--) for (int j = i - 1; j >= 0; j--) if (d[i] & (1LL << j)) d[i] ^= d[j]; for (int i = 0; i <= 60; i++) if (d[i]) p[cnt++] = d[i]; } long long kthquery(long long k) { int ret = 0; if (k >= (1LL << cnt)) return -1; for (int i = 60; i >= 0; i--) if (k & (1LL << i)) ret ^= p[i]; return ret; } 最终板子如下 struct L_B { long long d[61], p[61]; int cnt; L_B() { memset(d, 0, sizeof(d)); memset(p, 0, sizeof(p)); cnt = 0; } bool insert(long long val) { for (int i = 60; i >= 0; i--) if (val & (1LL << i)) { if (!d[i]) { d[i] = val; break; } val ^= d[i]; } return val > 0; } long long query_max() { long long ret = 0; for (int i = 60; i >= 0; i--) if ((ret ^ d[i]) > ret) ret ^= d[i]; return ret; } long long query_min() { for (int i = 0; i <= 60; i++) if (d[i]) return d[i]; return 0; } void rebuild() { for (int i = 60; i >= 0; i--) for (int j = i - 1; j >= 0; j--) if (d[i] & (1LL << j)) d[i] ^= d[j]; for (int i = 0; i <= 60; i++) if (d[i]) p[cnt++] = d[i]; } long long kthquery(long long k) { int ret = 0; if (k >= (1LL << cnt)) return -1; for (int i = 60; i >= 0; i--) if (k & (1LL << i)) ret ^= p[i]; return ret; } } L_B merge(const L_B & n1, const L_B & n2) { L_B ret = n1; for (int i = 60; i >= 0; i--) if (n2.d[i]) ret.insert(n1.d[i]); return ret; } BM杜教筛 #include <bits/stdc++.h> using namespace std; #define rep(i,a,n) for (long long i=a;i<n;i++) #define per(i,a,n) for (long long i=n-1;i>=a;i--) #define pb push_back #define mp make_pair #define all(x) (x).begin(),(x).end() #define fi first #define se second #define SZ(x) ((long long)(x).size()) typedef vector<long long> VI; typedef long long ll; typedef pair<long long,long long> PII; const ll mod=1e9+7; ll powmod(ll a,ll b) {ll res=1;a%=mod; assert(b>=0); for(;b;b>>=1){if(b&1)res=res*a%mod;a=a*a%mod;}return res;} // head long long _,n; namespace linear_seq { const long long N=10010; ll res[N],base[N],_c[N],_md[N]; vector<long long> Md; void mul(ll *a,ll *b,long long k) { rep(i,0,k+k) _c[i]=0; rep(i,0,k) if (a[i]) rep(j,0,k) _c[i+j]=(_c[i+j]+a[i]*b[j])%mod; for (long long i=k+k-1;i>=k;i--) if (_c[i]) rep(j,0,SZ(Md)) _c[i-k+Md[j]]=(_c[i-k+Md[j]]-_c[i]*_md[Md[j]])%mod; rep(i,0,k) a[i]=_c[i]; } long long solve(ll n,VI a,VI b) { // a 系数 b 初值 b[n+1]=a[0]*b[n]+... //printf("%d\n",SZ(b)); ll ans=0,pnt=0; long long k=SZ(a); assert(SZ(a)==SZ(b)); rep(i,0,k) _md[k-1-i]=-a[i];_md[k]=1; Md.clear(); rep(i,0,k) if (_md[i]!=0) Md.push_back(i); rep(i,0,k) res[i]=base[i]=0; res[0]=1; while ((1ll<<pnt)<=n) pnt++; for (long long p=pnt;p>=0;p--) { mul(res,res,k); if ((n>>p)&1) { for (long long i=k-1;i>=0;i--) res[i+1]=res[i];res[0]=0; rep(j,0,SZ(Md)) res[Md[j]]=(res[Md[j]]-res[k]*_md[Md[j]])%mod; } } rep(i,0,k) ans=(ans+res[i]*b[i])%mod; if (ans<0) ans+=mod; return ans; } VI BM(VI s) { VI C(1,1),B(1,1); long long L=0,m=1,b=1; rep(n,0,SZ(s)) { ll d=0; rep(i,0,L+1) d=(d+(ll)C[i]*s[n-i])%mod; if (d==0) ++m; else if (2*L<=n) { VI T=C; ll c=mod-d*powmod(b,mod-2)%mod; while (SZ(C)<SZ(B)+m) C.pb(0); rep(i,0,SZ(B)) C[i+m]=(C[i+m]+c*B[i])%mod; L=n+1-L; B=T; b=d; m=1; } else { ll c=mod-d*powmod(b,mod-2)%mod; while (SZ(C)<SZ(B)+m) C.pb(0); rep(i,0,SZ(B)) C[i+m]=(C[i+m]+c*B[i])%mod; ++m; } } return C; } long long gao(VI a,ll n) { VI c=BM(a); c.erase(c.begin()); rep(i,0,SZ(c)) c[i]=(mod-c[i])%mod; return solve(n,c,VI(a.begin(),a.begin()+SZ(c))); } }; int main() { while(~scanf("%lld", &n)) { printf("%lld\n",linear_seq::gao(VI{1,5,11,36,95,281,781,2245,6336,18061, 51205},n-1)); } }

模拟退火 HDU-2899 Strange Function 题目链接：[ HDU - 2899 ] 题目大意： 函数 F(x) = 6x7 + 8x6 + 7x3 + 5x2 - yx, 其中x的范围是0 ≤ x ≤ 100. 输入y值，输出F(x)的最小值 模拟退火算法 模拟退火就是类似于物体降温的概率，来进行多次搜索迭代 在迭代过程中，模拟退火算法随机选择下一个状态，有两种可能 新状态比原来状态更优，那么接受这个新状态 新状态更差，那么以一定的概率接受该状态，不过这个概率应该随着时间的推移逐渐降低 模拟退火算法的主要步骤如下： 设置一个初始的温度T 温度下降，状态转移。从当前温度按降温系数下降到下一个温度，在新的温度计算当前状态 如果温度降到设定的温度下界，程序停止 伪代码 eps = 1e-8; // 终止温度，接近0，用于控制精度 T = 100; // 初始温度，应该是高温，以100℃为例 delta = 0.98; // 降温系数，控制退火的快慢，小于1，以0.98为例 g(x); // 状态x时的评价函数，例如物理意义上的能量 now, next; // 当前状态和新状态 while (T > eps) { // 如果温度未降到eps g(next), g(now); // 计算能量 dE = g(next) - g(now); // 能量差 if (dE >= 0) // 新状态更优， 接受新状态 now = next; else if (exp(dE / T) > rand()) // 如果新状态更差，在一定概率下接受他，e^(dE/T) now = next; T *= delta; // 降温，退火过程模拟 } 常见问题 模拟退火算法在ACM中的典型问题有函数最值问题、TSP旅行商问题、最小圆覆盖、最小球覆盖等 题解代码 #include<bits/stdc++.h> using namespace std; const double eps = 1e-8; // 终止温度 double y; double func(double x) // 计算函数值 { return 6 * pow(x, 7.0) + 8 * pow(x, 6.0) + 7 * pow(x, 3.0) + 5 * pow(x, 2.0) - y * x; } double solve(void) { double T = 100; // 初始温度 double delta = 0.98; // 降温系数 double x = 50.0; // x的初始值 double now = func(x); // 计算初始函数值 double ans = now; // 返回值 while (T > eps) { int f[2] = {1, -1}; double newx = x + f[rand() % 2] * T; // 按概率改变x，随t的降温而减少 if (newx >= 0 && newx <= 100) { double next = func(newx); ans = min(ans, next); if (now - next > eps) { x = newx; now = next; } } T *= delta; } return ans; } int main(void) { int cas; scanf("%d", &cas); while (cas--) { scanf("%lf", &y); printf("%.4lf\n", solve()); } return 0; }

二维前缀和+差分 HDU6514 Monitor 题目链接：[ HDU - 6514 ] 题目大意：给你个n×m的区域，每个区域被选中标记为1，没被选中标记为0，然后给一些标记的区域，查询若干个区域问是否有0的 一道二维前缀和的题目，想通没什么难度 代码如下 #include <bits/stdc++.h> using namespace std; #define IO ios::sync_with_stdio(false), cin.tie(0), cout.tie(0) #define fin freopen("in.txt", "r", stdin) #define fout freopen("out.txt", "w", stdout) typedef long long ll; typedef unsigned long long ull; typedef pair<int, int> P; const int maxn = 1e7 + 10; const int INF = 0x3f3f3f3f; const double eps = 1e-5; int n, m, p, q; int a[maxn]; int getid(int i, int j) { if (i == 0 || j == 0 || i > n || j > m) return 0; return (i - 1) * m + j; } void add(int i, int j, int v) { int id = getid(i, j); if (!id) return; a[id] += v; } int main(void) { IO; while (~scanf("%d%d", &n, &m)) { memset(a, 0, sizeof(a)); int x1, y1, x2, y2; scanf("%d", &p); for (int i = 0; i < p; i++) { scanf("%d%d%d%d", &x1, &y1, &x2, &y2); add(x1, y1, 1); add(x2 + 1, y2 + 1, 1); add(x1, y2 + 1, -1); add(x2 + 1, y1, -1); } for (int i = 1; i <= n; i++) for (int j = 1; j <= m; j++) a[getid(i, j)] += a[getid(i - 1, j)] + a[getid(i, j - 1)] - a[getid(i - 1, j - 1)]; for (int i = 1; i <= n; i++) for (int j = 1; j <= m; j++) if (a[getid(i,j)] > 0) a[getid(i, j)] = 1; for (int i = 1; i <= n; i++) for (int j = 1; j <= m; j++) a[getid(i, j)] += a[getid(i - 1, j)] + a[getid(i, j - 1)] - a[getid(i - 1, j - 1)]; scanf("%d", &q); for (int i = 0; i < q; i++) { scanf("%d%d%d%d", &x1, &y1, &x2, &y2); int ans = a[getid(x2, y2)] - a[getid(x1 - 1, y2)] - a[getid(x2, y1 - 1)] + a[getid(x1 - 1, y1 - 1)]; if (ans == (x2 - x1 + 1) * (y2 - y1 + 1)) puts("YES"); else puts("NO"); } } return 0; }

树形DP POJ-2342 没有上司的舞会 题目链接：[POJ - 2342] 题目 某大学有N个职员，编号为1~N。他们之间有从属关系，也就是说他们的关系就像一棵以校长为根的树，父结点就是子结点的直接上司。现在有个周年庆宴会，宴会每邀请来一个职员都会增加一定的快乐指数Ri，但是呢，如果某个职员的上司来参加舞会了，那么这个职员就无论如何也不肯来参加舞会了。所以，请你编程计算，邀请哪些职员可以使快乐指数最大，求最大的快乐指数。 Input 第一行一个整数N。(1<=N<=6000) 接下来N行，第i+1行表示i号职员的快乐指数Ri。(-128<=Ri<=127) 接下来N行分别代表每个结点的权值范围从-128到127 接下来N-1行，每行输入一对整数L,K。表示K是L的直接上司。, 最后一行输入0 0 Output 对于每组测试数据，输出一个整数Ans，表示在不发生口角的情况下，乘务员最多可以清扫的垃圾数目。 Sample Input 5 2 1 36 9 80 69 85 Sample Output 201 解1 #include <bits/stdc++.h> using namespace std; int dp[6010][2]; int v[6010], a[6010]; vector<int> son[6010]; int n; void dfs(int x) { dp[x][0] = 0; dp[x][1] = a[x]; for (int i = 0; i < son[x].size(); i++) { int id = son[x][i]; dfs(id); dp[x][0] += max(dp[id][0], dp[id][1]); dp[x][1] += dp[id][0]; } } int main(void) { scanf("%d", &n); for (int i = 1; i <= n; i++) scanf("%d", &a[i]); int x, y; for (int i = 1; i < n; i++) { scanf("%d%d", &x, &y); son[y].push_back(x); v[x] = 1; } scanf("%d%d", &x, &y); int root = 0; for (int i = 1; i <= n; i++) { if (!v[i]) { root = i; break; } } dfs(root); cout << max(dp[root][1], dp[root][0]) << endl; getchar(), getchar(); return 0; } 解2 题解 这道题我们需要存储当前位置i开始往前m个位置的状态 在递推式dp[i][j] = max(dp[i - 1][j >> 1], dp[i - 1][(j >> 1) + (1 << m - 1)]) + w[i]中 dp[i - 1][j >> 1]、dp[i - 1][(j >> 1) + (1 << m - 1)]) + w[i]为再往前一个节点为0或者1的状态 #include<bits/stdc++.h> using namespace std; const int maxn = 1010; int dp[maxn][1050]; int w[maxn]; int n, m, q; int sum(int j) { int sum = 0; while (j) { if (j & 1) sum++; j >>= 1; } return sum; } int main(void) { memset(dp, 0, sizeof(dp)); cin >> n >> m >> q; for (int i = 1; i <= n; i++) scanf("%d", &w[i]); for (int i = 1; i <= n; i++) { for (int j = 0; j < 1 << m; j++) { if (sum(j) > q) continue; if (j & 1) dp[i][j] = max(dp[i - 1][j >> 1], dp[i - 1][(j >> 1) + (1 << m - 1)]) + w[i]; else dp[i][j] = max(dp[i - 1][j >> 1], dp[i - 1][(j >> 1) + (1 << m - 1)]); } } int sum = -1; for (int i = 0; i < 1 << m; i++) sum = max(sum, dp[n][i]); cout << sum << endl; return 0; }

状压DP Hiho-1044 状态压缩 题目链接：[HihoCoder - 1044] 题目 小Hi和小Ho在兑换到了喜欢的奖品之后，便继续起了他们的美国之行，思来想去，他们决定乘坐火车前往下一座城市——那座城市即将举行美食节！ 但是不幸的是，小Hi和小Ho并没有能够买到很好的火车票——他们只能够乘坐最为破旧的火车进行他们的旅程。 不仅如此，因为美食节的吸引，许多人纷纷踏上了和小Hi小Ho一样的旅程，于是有相当多的人遭遇到了和小Hi小Ho一样的情况——这导致这辆车上的人非常非常的多，以至于都没有足够的位置能让每一个人都有地方坐下来。 小Hi和小Ho本着礼让他们的心情——当然还因为本来他们买的就是站票，老老实实的呆在两节车厢的结合处。他们本以为就能够这样安稳抵达目的地，但事与愿违，他们这节车厢的乘务员是一个强迫症，每隔一小会总是要清扫一次卫生，而时值深夜，大家都早已入睡，这种行为总是会惊醒一些人。而一旦相邻的一些乘客被惊醒了大多数的话，就会同乘务员吵起来，弄得大家都睡不好。 将这一切看在眼里的小Hi与小Ho决定利用他们的算法知识，来帮助这个有着强迫症的乘务员——在不与乘客吵起来的前提下尽可能多的清扫垃圾。 小Hi和小Ho所处的车厢可以被抽象成连成一列的N个位置，按顺序分别编号为1..N，每个位置上都有且仅有一名乘客在休息。同时每个位置上都有一些垃圾需要被清理，其中第i个位置的垃圾数量为Wi。乘务员可以选择其中一些位置进行清理，但是值得注意的是，一旦有编号连续的M个位置中有超过Q个的位置都在这一次清理中被选中的话（即这M个位置上的乘客有至少Q+1个被惊醒了），就会发生令人不愉快的口角。而小Hi和小Ho的任务是，计算选择哪些位置进行清理，在不发生口角的情况下，清扫尽可能多的垃圾。 Input 每个测试点（输入文件）有且仅有一组测试数据。 每组测试数据的第一行为三个正整数N、M和Q，意义如前文所述。 每组测试数据的第二行为N个整数，分别为W1到WN，代表每一个位置上的垃圾数目。 对于100%的数据，满足N<=1000, 2<=M<=10,1<=Q<=M, Wi<=100 Output 对于每组测试数据，输出一个整数Ans，表示在不发生口角的情况下，乘务员最多可以清扫的垃圾数目。 Sample Input 5 2 1 36 9 80 69 85 Sample Output 201 题解 这道题我们需要存储当前位置i开始往前m个位置的状态 在递推式dp[i][j] = max(dp[i - 1][j >> 1], dp[i - 1][(j >> 1) + (1 << m - 1)]) + w[i]中 dp[i - 1][j >> 1]、dp[i - 1][(j >> 1) + (1 << m - 1)]) + w[i]为再往前一个节点为0或者1的状态 #include<bits/stdc++.h> using namespace std; const int maxn = 1010; int dp[maxn][1050]; int w[maxn]; int n, m, q; int sum(int j) { int sum = 0; while (j) { if (j & 1) sum++; j >>= 1; } return sum; } int main(void) { memset(dp, 0, sizeof(dp)); cin >> n >> m >> q; for (int i = 1; i <= n; i++) scanf("%d", &w[i]); for (int i = 1; i <= n; i++) { for (int j = 0; j < 1 << m; j++) { if (sum(j) > q) continue; if (j & 1) dp[i][j] = max(dp[i - 1][j >> 1], dp[i - 1][(j >> 1) + (1 << m - 1)]) + w[i]; else dp[i][j] = max(dp[i - 1][j >> 1], dp[i - 1][(j >> 1) + (1 << m - 1)]); } } int sum = -1; for (int i = 0; i < 1 << m; i++) sum = max(sum, dp[n][i]); cout << sum << endl; return 0; }

数位DP 回文序列 POJ-3280 Cheapest Palindrome 题目链接：[ POJ - 3280 ] 题目大意 给定字符串s，长度为m，由n个小写字母组成。在s的任意位置增删字母，把它变成回文串，增删特定字母的花费不同，求最小花费 思路 定义状态 dp[i][j] 表示字符串s的子区间 s[i, j] 变成回文的最小花费 那么每次有三种情况: 如果 s[i] == s[j], 那么 dp[i][j] = dp[i + 1][j - 1]。 如果 dp[i + 1][j] 是回文串, 那么 dp[i][j] = dp[i + 1][j] + w[i]。 如果 dp[i][j - 1] 是回文串, 那么 dp[i][j] = dp[i][j - 1] + w[j]。 2, 3情况的状态转移方程就是 dp[i][j] = min(dp[i + 1][j] + w[i], dp[i][j - 1] + w[j]) 题解代码 #include<bits/stdc++.h> using namespace std; int w[30], dp[2010][2010]; char s[2010], ch; int n, m; int main(void) { int x, y; while (cin >> n >> m) { cin >> s; for (int i = 0; i < n; i++) { cin >> ch >> x >> y; w[ch - 'a'] = min(x, y); } for (int i = m - 1; i >= 0; i--) { // i是子区间的起点 for (int j = i + 1; j < m; j++) { // j是子区间的终点 if (s[i] == s[j]) dp[i][j] = dp[i + 1][j - 1]; else dp[i][j] = min(dp[i + 1][j] + w[s[i] - 'a'], dp[i][j - 1] + w[s[j] - 'a']); } } cout << dp[0][m - 1] << endl; } return 0; }

数位DP 石子合并 模板题 题目 有n堆石子排成一排，每堆石子有一定的数量，将n堆石子合并成一堆。合并的规则是每次只能合并相邻的两堆石子，合并的花费为这两堆石子的总数，求最小花费。 题解代码 #include<bits/stdc++.h> using namespace std; int sum[1010]; int dp[1010][1010], s[1010][1010]; int n, x; int main(void) { while (~scanf("%d", &n)) { sum[0] = 0; for (int i = 1; i <= n; i++) { scanf("%d", &x); sum[i] = sum[i - 1] + x; } for (int i = 1; i <= n; i++) { dp[i][i] = 0; s[i][i] = i; } for (int len = 1; len < n; len++) { for (int i = 1; i <= n - len; i++) { int j = i + len; dp[i][j] = INF; for (int k = s[i][j - 1]; k <= s[i + 1][j]; k++) { if (dp[i][k] + dp[k + 1][j] + sum[j] - sum[i -1] < dp[i][j]) { dp[i][j] = dp[i][k] + dp[k + 1][j] + sum[j] - sum[i - 1]; s[i][j] = k; } } } } cout << dp[1][n] << endl; } return 0; }

线段树+扫描线 HDU-3642 Get The Treasury 题目链接：[ HDU - 3642 ] 题目大意：给你n个立方体的左下角坐标和右上角坐标，求立方体相交至少相交的面积。 题解 将Z轴离散化，把面积当作二维的底和高 代码如下 #include<bits/stdc++.h> using namespace std; const int MAXN=2010; struct Node { int l,r; int lf,rf;//实际的左右端点 int c; int once,twice,more; }segTree[MAXN*3]; int y[MAXN]; int z[MAXN]; struct Line { int x; int y1,y2; int z1,z2;//这两个是枚举的时候判断使用的 int f; }line[MAXN]; bool cmp(Line a,Line b) { return a.x<b.x; } void Build(int i,int l,int r) { segTree[i].l=l; segTree[i].r=r; segTree[i].lf=y[l]; segTree[i].rf=y[r]; segTree[i].c=0; segTree[i].once=segTree[i].twice=segTree[i].more=0; if(r==l+1)return; int mid=(l+r)>>1; Build(i<<1,l,mid); Build((i<<1)|1,mid,r); } void push_up(int i) { if(segTree[i].c>2) { segTree[i].more=segTree[i].rf-segTree[i].lf; segTree[i].once=segTree[i].twice=0; } else if(segTree[i].c==2) { if(segTree[i].l+1==segTree[i].r)//叶子结点 { segTree[i].more=0; segTree[i].twice=segTree[i].rf-segTree[i].lf; segTree[i].once=0; return; } segTree[i].more=segTree[i<<1].once+segTree[i<<1].twice+segTree[i<<1].more +segTree[(i<<1)|1].once+segTree[(i<<1)|1].twice+segTree[(i<<1)|1].more; segTree[i].twice=segTree[i].rf-segTree[i].lf-segTree[i].more; segTree[i].once=0; } else if(segTree[i].c==1) { if(segTree[i].l+1==segTree[i].r) { segTree[i].more=0; segTree[i].twice=0; segTree[i].once=segTree[i].rf-segTree[i].lf; return; } segTree[i].more=segTree[i<<1].more+segTree[i<<1].twice +segTree[(i<<1)|1].more+segTree[(i<<1)|1].twice; segTree[i].twice=segTree[i<<1].once+segTree[(i<<1)|1].once; segTree[i].once=segTree[i].rf-segTree[i].lf-segTree[i].more-segTree[i].twice; } else { if(segTree[i].l+1==segTree[i].r) { segTree[i].more=segTree[i].once=segTree[i].twice=0; return; } segTree[i].more=segTree[i<<1].more+segTree[(i<<1)|1].more; segTree[i].twice=segTree[i<<1].twice+segTree[(i<<1)|1].twice; segTree[i].once=segTree[i<<1].once+segTree[(i<<1)|1].once; } } void update(int i,Line e) { if(segTree[i].lf>=e.y1 && segTree[i].rf<=e.y2) { segTree[i].c+=e.f; push_up(i); return; } if(e.y2<=segTree[i<<1].rf) update(i<<1,e); else if(e.y1>=segTree[(i<<1)|1].lf) update((i<<1)|1,e); else { Line temp=e; temp.y2=segTree[i<<1].rf; update(i<<1,temp); temp=e; temp.y1=segTree[(i<<1)|1].lf; update((i<<1)|1,temp); } push_up(i); } Line temp[MAXN]; int main(void) { int T; int n; int x1,y1,z1,x2,y2,z2; scanf("%d",&T); int iCase=0; while(T--) { iCase++; scanf("%d",&n); int t=0; for(int i=0;i<n;i++) { scanf("%d%d%d%d%d%d", &x1, &y1, &z1, &x2, &y2, &z2); line[t].x=x1; line[t].y1=y1; line[t].y2=y2; line[t].z1=z1; line[t].z2=z2; line[t].f=1; y[t]=y1; z[t++]=z1; line[t].x=x2; line[t].y1=y1; line[t].y2=y2; line[t].z1=z1; line[t].z2=z2; line[t].f=-1; y[t]=y2; z[t++]=z2; } sort(line,line+t,cmp); sort(y,y+t); int t1=unique(y,y+t)-y; Build(1,0,t1-1); sort(z,z+t); int t2=unique(z,z+t)-z; long long ans=0; long long area=0; for(int i=0;i<t2-1;i++) { int m=0; for(int j=0;j<t;j++) if(line[j].z1<=z[i]&&line[j].z2>z[i]) temp[m++]=line[j]; area=0; update(1,temp[0]); for(int j=1;j<m;j++) { area+=(long long)segTree[1].more*(temp[j].x-temp[j-1].x); update(1,temp[j]); } ans+=area*(z[i+1]-z[i]); } printf("Case %d: %I64d\n",iCase,ans); } return 0; }

线段树+扫描线 HDU-1542 Picture 题目链接：[ HDU - 1828 ] 题目大意：给你n个矩形的左下角坐标和右上角坐标，求外周长。 题解 扫描线 + 离散化板子 将横竖两次扫描简化 代码如下 #include<bits/stdc++.h> using namespace std; const int maxn = 10010; const int maxN = 80010; const int INF = 0x3f3f3f3f; inline int read(void) { int w = 1, d = 0; char ch = getchar(); while (ch < '0' || ch > '9') { if (ch == '-') w = -1; ch = getchar(); } while (ch >= '0' && ch <= '9') { d = d * 10 + ch - '0'; ch = getchar(); } return w * d; } struct Node { int l, r, cnt; int lc, rc; int len, num; } tree[maxN]; struct node { int x1, x2, y, flag; void init(int l, int r, int h, int f) { x1 = l, x2 = r, y = h, flag = f; } bool operator < (const node& a) { return y < a.y; } } line[maxn]; int n, x1, y11, x2, y2; void build_tree(int id, int l, int r) { tree[id].l = l; tree[id].r = r; tree[id].cnt = 0; tree[id].len = 0; tree[id].lc = 0; tree[id].rc = 0; tree[id].num = 0; if (l == r) return; int mid = (l + r) >> 1; build_tree(id << 1, l, mid); build_tree(id << 1 | 1, mid + 1, r); } void getlen(int i) { if (tree[i].cnt) { tree[i].len = tree[i].r - tree[i].l + 1; tree[i].lc = tree[i].rc = 1; tree[i].num = 1; }else if (tree[i].l == tree[i].r) { tree[i].len = 0; tree[i].lc = tree[i].rc = 0; tree[i].num = 0; } else { tree[i].len = tree[i << 1].len + tree[i << 1 | 1].len; tree[i].lc = tree[i << 1].lc; tree[i].rc = tree[i << 1 | 1].rc; tree[i].num = tree[i << 1].num + tree[i << 1 | 1].num - (tree[i << 1].rc & tree[i << 1 | 1].lc); } } void update(int id, int l, int r, int v) { if (tree[id].l == l && tree[id].r == r) { tree[id].cnt += v; getlen(id); return; } int mid = (tree[id].l + tree[id].r) >> 1; if (r <= mid) update(id << 1, l, r, v); else if (l > mid) update(id << 1 | 1, l, r, v); else { update(id << 1, l, mid, v); update(id << 1 | 1, mid + 1, r, v); } getlen(id); } int main(void) { ios::sync_with_stdio(false); while (~scanf("%d", &n)) { int minn = INF, maxx = -INF; for (int i = 0; i < n; i++) { x1 = read(), y11 = read(), x2 = read(), y2 = read(); maxx = max(maxx, max(x1, x2)); minn = min(minn, min(x1, x2)); line[i].init(x1, x2, y11, 1); line[i + n].init(x1, x2, y2, -1); } sort(line, line + 2 * n); build_tree(1, minn, maxx - 1); long long ans = 0, last = 0; for (int i = 0; i < 2 * n; i++) { update(1, line[i].x1, line[i].x2 - 1, line[i].flag); ans += abs(tree[1].len - last); ans += (line[i + 1].y - line[i].y) * tree[1].num * 2; last = tree[1].len; } cout << ans << endl; } return 0; }

线段树+扫描线 HDU-1542 Atlantis 题目链接：[ HDU - 1542 ] 题目大意：给你n个矩形的左下角坐标和右上角坐标，求矩形相交的面积。 题解 扫描线 + 离散化板子 代码如下 #include<cstdio> #include<cstring> #include<algorithm> using namespace std; double dif_x[500]; struct node { double x1,x2,y; int flag; void init(double l, double r, double h, int key) { x1 = l; x2 = r; y = h; flag = key; } }line[500]; bool cmp(node a, node b) { return a.y < b.y; } struct node2 { int l, r, cnt; double len; }tree[1000]; void build_tree(int id, int l, int r) { tree[id].l = l; tree[id].r = r; tree[id].cnt = 0; tree[id].len = 0; if(l == r) return; int mid = (r + l) >> 1; build_tree(id << 1, l, mid); build_tree((id << 1) + 1, mid + 1, r); } void getlen(int id) { if(tree[id].cnt >= 1) tree[id].len = dif_x[tree[id].r+1] - dif_x[tree[id].l]; else tree[id].len = tree[id<<1].len + tree[(id<<1)|1].len; } void update(int id,int l,int r,int v) { if(tree[id].l==l && tree[id].r==r) { tree[id].cnt += v; getlen(id); return; } int mid = (tree[id].l + tree[id].r) >> 1; if(r <= mid) update(id << 1,l,r,v); else if(l > mid) update((id << 1) + 1, l, r, v); else { update(id << 1, l, mid, v); update((id << 1) + 1, mid + 1, r, v); } getlen(id); } int mySearch(double p, int l, int r) { while(l <= r) { int mid = (l + r) >> 1; if(dif_x[mid] == p) return mid; if(dif_x[mid] < p) l = mid + 1; else r = mid - 1; } } int main(void) { int n, noc = 0; while (~scanf("%d", &n)) { if (n == 0) break; noc++; double x1, y1, x2, y2; int line_num = 0; for(int i = 0; i < n; i++) { scanf("%lf%lf%lf%lf", &x1, &y1, &x2, &y2); line[line_num].init(x1, x2, y1, 1); dif_x[line_num++] = x1; line[line_num].init(x1, x2, y2, -1); dif_x[line_num++] = x2; } sort(line, line + line_num, cmp); sort(dif_x, dif_x + line_num); int dif_x_num = unique(dif_x, dif_x + line_num) - dif_x; build_tree(1, 0, dif_x_num - 1); double ans = 0; for(int i = 0; i < line_num - 1; i++) { int line_l = mySearch(line[i].x1, 0, dif_x_num - 1); int line_r = mySearch(line[i].x2, 0, dif_x_num-1) - 1; update(1, line_l, line_r, line[i].flag); ans += tree[1].len * (line[i+1].y - line[i].y); } printf("Test case #%d\n", noc); printf("Total explored area: %.2lf\n\n", ans); } return 0; }

线段树+扫描线 HDU-1542 覆盖的面积 题目链接：[ HDU - 1255 ] 题目大意：给你n个矩形的左下角坐标和右上角坐标，求矩形相交至少覆盖两次以上的面积。 题解 代码如下 #include<bits/stdc++.h> using namespace std; const int maxn = 2010; double dif_x[maxn]; //记录不同的x坐标 int n, t = 0; struct node { double x1, x2, y; int flag; void init(double l, double r, double h, int key) { x1 = l, x2 = r, y = h, flag = key; } bool operator < (const node& a) const { return y < a.y; } }line[maxn]; struct Node { int l, r, cnt; double len; double len2; }tree[maxn * 4]; void build_tree(int id, int l, int r) { tree[id].l = l; tree[id].r = r; tree[id].cnt = 0; tree[id].len = 0; tree[id].len2 = 0; if (l == r) return; int mid = (r + l) >> 1; build_tree(id << 1, l, mid); build_tree(id << 1 | 1, mid + 1, r); } void getlen(int id) { if (tree[id].cnt) //如果该段被覆盖那么就直接由dif_x数组获得长度 tree[id].len = dif_x[tree[id].r + 1] - dif_x[tree[id].l]; else if (tree[id].l == tree[id].r) tree[id].len = 0; else //如果没有被覆盖，那么应该是由左右孩子的和 tree[id].len = tree[id << 1].len + tree[id << 1 | 1].len; if (tree[id].cnt >= 2) tree[id].len2 = dif_x[tree[id].r + 1] - dif_x[tree[id].l]; else if (tree[id].l == tree[id].r) tree[id].len2 = 0; else if (tree[id].cnt == 1) tree[id].len2 = tree[id << 1].len + tree[id << 1 | 1].len; else tree[id].len2 = tree[id << 1].len2 + tree[id << 1 | 1].len2; } void update(int id, int l, int r, int v) { if (tree[id].l == l && tree[id].r == r) { tree[id].cnt += v; getlen(id); return; } int mid = (tree[id].l + tree[id].r) >> 1; if (r <= mid) update(id << 1, l, r, v); else if (l > mid) update(id << 1 | 1, l, r, v); else { update(id << 1, l, mid, v); update(id << 1 | 1, mid + 1, r, v); } getlen(id); } int main(void) { ios::sync_with_stdio(false); scanf("%d", &t); while (t--) { scanf("%d", &n); int line_num = 0; //一共有多少条横线 double x1, y11, x2, y2; for (int i = 0; i < n; i++) { scanf("%lf%lf%lf%lf", &x1, &y11, &x2, &y2); line[line_num].init(x1, x2, y11, 1); dif_x[line_num++] = x1; line[line_num].init(x1, x2, y2, -1); dif_x[line_num++] = x2; } sort(line, line + line_num); sort(dif_x, dif_x + line_num); //对dif_x去重，要先排个序，这样更方便 int dif_x_num = unique(dif_x, dif_x + line_num) - dif_x; //dif_x_num表示去重后不同的x坐标的数量 build_tree(1, 0, dif_x_num - 1); double ans = 0.00; for (int i = 0; i < line_num - 1; i++) { int ll = lower_bound(dif_x, dif_x + dif_x_num, line[i].x1) - dif_x; int rr = lower_bound(dif_x, dif_x + dif_x_num, line[i].x2) - dif_x - 1; update(1, ll, rr, line[i].flag); ans += tree[1].len2 * (line[i + 1].y - line[i].y); } printf("%.2lf\n", ans + 0.000001); } return 0; }

树状数组 POJ-2155 Matrix 题目链接：[ POJ - 2155 ] 题目大意 对一个n∗n的矩阵： 格式C x1 y1 x2 y2，表示将左上角为(x1,y1)，右下角为(x2,y2)的矩阵全部取反，即0变1，1变0. Q x y，表示查询位置(x,y)的值. 设询问次数为t，则10^3, 1 ≤ t ≤ 10^5, 数据组数 ≤ 10. 题解 树状数组模板题，此题建立完整的二维树状数组板子 代码如下 #include<iostream> #include<algorithm> #include<cstring> #include<cstdio> #include<string> #include<vector> using namespace std; typedef long long ll; ll a[1010][1010]; ll t, m, n; inline ll lowbit(ll x) { return x & -x; } ll sum(ll x, ll y) { ll res = 0; for (ll i = x; i > 0; i -= lowbit(i)) for (ll j = y; j > 0; j -= lowbit(j)) res += a[i][j]; return res; } void add(ll x, ll y, ll v) { for (ll i = x; i <= n; i += lowbit(i)) for (ll j = y; j <= n; j += lowbit(j)) a[i][j] += v; } int main(void) { cin >> t; while (t--) { memset(a, 0, sizeof(a)); scanf("%lld%lld", &n, &m); getchar(); char ch; ll x1, y1, x2, y2; while (m--) { scanf("%c %lld %lld", &ch, &x1, &y1); if (ch == 'Q') printf("%lld\n", sum(x1, y1) % 2); if (ch == 'C') { scanf("%lld%lld", &x2, &y2); add(x1, y1, 1); add(x2 + 1, y1, 1); add(x1, y2 + 1, 1); add(x2 + 1, y2 + 1, 1); } getchar(); } if (t) cout << endl; } return 0; }

树状数组 POJ-3468 A Simple Problem with Integers 题目链接：[ POJ - 3468 ] 题目大意 给一个长度为n的数列，有Q次操作Q代表查询区间 a b之间的累加和，操作C代表将a-b区间的所有数加上c 题解 树状数组模板题，此题建立完整的一维树状数组板子 代码如下 #include<iostream> #include<algorithm> #include<cstdio> #include<string> #include<vector> using namespace std; typedef long long ll; ll n, m, k; inline ll lowbit(ll x) { return x & -x; } struct BIT { ll a[500050]; void add(ll x, ll v) { for (x; x > 0; x -= lowbit(x)) a[x] += v; } void addinv(ll x, ll v) { for (ll i = x; i <= n; i += lowbit(i)) a[i] += x * v; } ll Sum(ll x) { ll res = 0; for(x; x <= n; x += lowbit(x)) res += a[x]; return res; } ll Suminv(ll x) { ll res = 0; for (x; x > 0; x -= lowbit(x)) res += a[x]; return res; } }bit, bitinv; ll sum(ll x) { if (x) return bit.Sum(x) * x + bitinv.Suminv(x - 1); else return 0; } void add(ll l, ll r, ll c) { bit.add(r, c); bitinv.addinv(r, c); if (l > 1) { bit.add(l - 1, -c); bitinv.addinv(l - 1, -c); } } int main(void) { scanf("%lld%lld", &n, &m); for (ll i = 1; i <= n; i++) { scanf("%lld", &k); add(i, i, k); } getchar(); char ch; ll l, r, c; while (m--) { scanf("%c %lld %lld", &ch, &l, &r); if (ch == 'Q') printf("%lld\n", sum(r) - sum(l - 1)); if (ch == 'C') { scanf("%lld", &c); add(l, r, c); } getchar(); } return 0; }