zoukankan      html  css  js  c++  java
  • 二维几何基础

    1. 点、线、凸边形

      1 /*******************************************************
      2                     二维几何基础
      3 【注意】数组下标从1开始。
      4 *******************************************************/
      5 
      6 #include <iostream>
      7 #include <cstdio>
      8 #include <cstdlib>
      9 #include <cmath>
     10 #include <iomanip>
     11 #include <algorithm>
     12 #include <string>
     13 #define re register
     14 #define il inline
     15 #define ll long long
     16 #define ld long double
     17 using namespace std;
     18 const ll MAXN = 1e5+5;
     19 const ll INF = 1e9;
     20 const ld EPS = 1e-9;
     21 
     22 //点坐标
     23 struct POINT
     24 {
     25     ld x, y;
     26     POINT() : x(0), y(0) {}
     27     POINT(ld _x, ld _y) : x(_x), y(_y) {}
     28 };
     29 //向量
     30 typedef POINT VECTOR;
     31 
     32 POINT xy[MAXN];     //顺时针多边形顶点存储
     33 POINT xyB[MAXN];    //判断点存储
     34 
     35 
     36 //符号函数
     37 ll sgn(ld x)    {return x < -EPS ? -1 : x > EPS;}
     38 //向量+向量=向量,点+向量=点
     39 VECTOR operator + (VECTOR a, VECTOR b)  {return {a.x+b.x,a.y+b.y};}
     40 //向量-向量=向量,点-点=向量
     41 VECTOR operator - (VECTOR a, VECTOR b)  {return {a.x-b.x,a.y-b.y};}
     42 //向量*数=向量
     43 VECTOR operator * (VECTOR a, ld k)  {return {a.x*k,a.y*k};}
     44 VECTOR operator * (ld k, VECTOR a)  {return {a.x*k,a.y*k};}
     45 //向量/数=向量
     46 VECTOR operator / (VECTOR a, ld k)  {return {a.x/k,a.y/k};}
     47 //向量==向量,点==点
     48 bool operator == (const VECTOR a, const VECTOR b)   {return !sgn(a.x-b.x) && !sgn(a.y-b.y);}
     49 //向量偏序关系(先x轴再y轴)
     50 bool operator < (const VECTOR a, const VECTOR b)    {return !sgn(a.x-b.x) ? a.y < b.y : a.x < b.x;}
     51 //向量旋转(逆时针)
     52 VECTOR rot(VECTOR a, ld sita)   {return {a.x*cos(sita)-a.y*sin(sita),a.x*sin(sita)+a.y*cos(sita)};}
     53 //取下端点
     54 POINT min(POINT p1, POINT p2)   {return p1.y<p2.y ? p1:p2;}
     55 //取上端点
     56 POINT max(POINT p1, POINT p2)   {return p1.y<p2.y ? p2:p1;}
     57 //点积
     58 ld dot(POINT p1, POINT p2, POINT p) {return (p.x-p1.x)*(p2.x-p1.x)+(p.y-p1.y)*(p2.y-p1.y);} //(三点式:共p1端点)
     59 ld dot(VECTOR a, VECTOR b)      {return a.x*b.x+a.y*b.y;}   //向量式
     60 //叉积
     61 ld cross(POINT p1, POINT p2, POINT p)   {return (p.x-p1.x)*(p2.y-p1.y)-(p.y-p1.y)*(p2.x-p1.x);} //(三点式:共p1端点)
     62 ld cross(VECTOR a, VECTOR b)    {return a.x*b.y-a.y*b.x;}   //向量式(aXb)
     63 //向量长度
     64 ld vlen(VECTOR a) {return sqrt(dot(a,a));}
     65 //向量夹角余弦值
     66 ld vcos(VECTOR a, VECTOR b) {return dot(a,b)/(vlen(a)*vlen(b));}
     67 //向量夹角正弦值
     68 ld vsin(VECTOR a, VECTOR b) {return fabs(cross(a,b))/(vlen(a)*vlen(b));}
     69 //向量极角
     70 ld vagl(VECTOR a, VECTOR b) {return acos(dot(a,b)/(vlen(a)*vlen(b)));}  //向量间
     71 ld vagl(VECTOR a) {return acos(a.x/vlen(a));}
     72 //求直线斜率【注意】确保斜率存在
     73 ld slope(VECTOR a)  {return a.y/a.x;}   //向量式
     74 ld slope(POINT p, POINT q)  {return (p.y-q.y)/(p.x-q.x);}   //两点式
     75 //单位向量【注意】确保非零向量
     76 VECTOR unit(VECTOR a)   {return a/vlen(a);}
     77 //两直线交点
     78 POINT intersectline(POINT p, VECTOR v, POINT q, VECTOR w)   {return p+v*cross(w,p-q)/cross(v,w);}   //(参数式:P=P0+t*v,P0为直线上某一点,v为直线的方向向量)
     79 //点在直线上的投影
     80 POINT proline(POINT a, POINT b, POINT p)    {return a+(b-a)*(dot(b-a,p-a)/dot(b-a,b-a));}
     81 //点关于直线的对称点
     82 POINT refline(POINT a, POINT b, POINT p)    {return proline(a,b,p)*2-p;}
     83 //判断两直线是否平行
     84 bool parallel(POINT p1, POINT p2, POINT q1, POINT q2)   {return !sgn(cross(p2-p1,q2-q1)) && sgn(cross(p1-q1,p2-q1));}
     85 //判断两直线重合
     86 bool superposition(POINT p1, POINT p2, POINT q1, POINT q2)  {return !sgn(cross(p2-p1,q2-q1)) && !sgn(cross(p1-q1,p2-q1));}
     87 //点到直线距离
     88 ld disline(POINT a, POINT b, POINT p)   {return fabs(cross(b-a,p-a))/vlen(b-a);}    //不取绝对值得到的是有向距离
     89 //点到线段距离(两种情况:点的投影在线段上,则为垂直距离;点的投影不在线段上,则为到两端距离的最小值)
     90 ld disseg(POINT a, POINT b, POINT p)    {return a==b ? vlen(p-a):dot(b-a,p-a)<0 ? vlen(p-a):dot(b-a,p-b)>0 ? vlen(p-b):fabs(cross(b-a,p-a))/vlen(b-a);}
     91 //线段相交判断(严格)
     92 bool intersectseg(POINT p1, POINT p2, POINT q1, POINT q2)   {return cross(p1-q1,q2-q1)*cross(p2-q1,q2-q1)<0 && cross(q1-p1,p2-p1)*cross(q2-p1,p2-p1)<0;}
     93 //判断点p(x,y)是否在线段p1p2上(两种情况:1.包括端点;2.不包含端点)
     94 bool onseg(POINT p1, POINT p2, POINT p) {return !sgn(cross(p1,p2,p)) && sgn(dot(p,p1,p2))<=0;}  //包含端点
     95 bool onseg_strict(POINT p1, POINT p2, POINT p)  {return !sgn(cross(p1,p2,p)) && sgn(dot(p,p1,p2))<0;}   //不包含端点
     96 
     97 //点射法判断点是否在多边形内部(边界也算在内部)
     98 //复杂度O(n)(不论凹凸)
     99 bool inpolygon_dot(ld x, ld y, ll n)
    100 {
    101     ll cnt = 0;
    102     for(re ll i = 0; i < n; ++i)
    103     {
    104         POINT p1 = min(xy[i+1], xy[(i+1)%n+1]);    //取下端点
    105         POINT p2 = max(xy[i+1], xy[(i+1)%n+1]);    //取上端点
    106         //特判点在线段上
    107         if(onseg(p1,p2,{x,y}))    return true;
    108         //从点(x,y)向x反方向引出射线
    109         //计算点射到的多边形的边数边数
    110         if(sgn(p1.y-y)<0 && sgn(y-p2.y) <=0 && sgn(cross(p1,p2,{x,y}))>0)   ++cnt;
    111     }
    112     if(cnt%2)   return true;
    113     else    return false;
    114 }
    115 
    116 //二分法判断点是否在多边形内部(不包括边界)
    117 //复杂度O(logn)(要求凸多边形)
    118 bool inpolygon_bis(ld x, ld y, ll n)
    119 {
    120     POINT p = {x,y};
    121     ll l = 2, r = n;
    122     //判断是否在幅角最小和最大的两条边上
    123     if(onseg(xy[1],xy[l],p) || onseg(xy[1],xy[n],p))    return false;
    124     //二分法判断是否在多边形内部
    125     while(l < r)
    126     {
    127         ll mid = (l+r)>>1;  //【注意】如此取中点最终:r==l或r==l-1(只有此两种情况)
    128         ll d = sgn(cross(xy[mid]-xy[1],p-xy[1]));
    129         if(d < 0)   l = mid+1;
    130         else if(d > 0)  r = mid-1;
    131         else
    132         {
    133             if(onseg_strict(xy[1],xy[mid],p))
    134             {
    135                 return true;
    136             }
    137             else
    138             {
    139                 return false;
    140             }
    141         }
    142     }
    143     //判断在最终二分得到的对角线两端的边是否都满足严格在内的条件
    144     if(l >= r && (sgn(cross(xy[l]-xy[l-1],p-xy[l-1]))>=0 || sgn(cross(xy[l%n+1]-xy[l],p-xy[l]))>=0))
    145     {
    146         return false;
    147     }
    148     return true;
    149 }
    150 
    151 //计算多边形面积(凹凸均可)
    152 ld polygonarea(ll n)
    153 {
    154     ld aera = 0;
    155     for(re ll i = 0; i < n; ++i)
    156     {
    157         aera += cross(xy[i+1], xy[(i+1)%n+1]);
    158     }
    159     //计算出来的aera为有向面积
    160     return fabs(aera)/2;
    161 }
    162 
    163 int main()
    164 {
    165     std::ios::sync_with_stdio(false);
    166     //判断直线平行/重合/相交
    167     ll n;
    168     cin >> n;
    169     cout << "INTERSECTING LINES OUTPUT" << endl;
    170     for(re ll i = 1; i <= n; ++i)
    171     {
    172         POINT p1, p2, p3, p4;
    173         cin >> p1.x >> p1.y >> p2.x >> p2.y;
    174         cin >> p3.x >> p3.y >> p4.x >> p4.y;
    175         if(parallel(p1,p2,p3,p4))   cout << "NONE" << endl;
    176         else if(superposition(p1,p2,p3,p4)) cout << "LINE" << endl;
    177         else
    178         {
    179             POINT q = intersectline(p1,p2-p1,p3,p4-p3);
    180             cout << "POINT " << fixed << setprecision(2) << q.x << " " << q.y << endl;
    181         }
    182     }
    183     cout << "END OF OUTPUT" << endl;
    184     return 0;
    185     /*
    186     //计算多边形的面积
    187     ll n;
    188     while(true)
    189     {
    190         std::cin >> n;
    191         if(!n)  break;
    192         for(re ll i = 1; i <= n; ++i)
    193         {
    194             std::cin >> xy[i].x;
    195             std::cin >> xy[i].y;
    196         }
    197         std::cout << std::fixed << std::setprecision(1) << polygonarea(n) << std::endl;
    198     }
    199     return 0;
    200     */
    201     /*
    202     //判断点是否在多边形内(不包括边界)
    203     ll n;
    204     std::cin >> n;
    205     for(re ll i = 1; i <= n; ++i)
    206     {
    207         std::cin >> xy[i].x;
    208         std::cin >> xy[i].y;
    209     }
    210     ll m;
    211     cin >> m;
    212     for(re ll i = 1; i <= m; ++i)
    213     {
    214         cin >> xyB[i].x;
    215         cin >> xyB[i].y;
    216     }
    217     for(re ll i = 1; i <= m; ++i)
    218     {
    219         if(!inpolygon_bis(xyB[i].x, xyB[i].y, n))
    220         {
    221             printf("NO
    ");
    222             return 0;
    223         }
    224     }
    225     printf("YES
    ");
    226     return 0;
    227     */
    228     /*
    229     //判断点是否在多边形内(包括边界)
    230     ll N, M, k = 0;
    231     while(!(std::cin >> N).eof())
    232     {
    233         if(!N)  break;
    234         std::cin >> M;
    235         if(k)   printf("
    ");
    236         ++k;
    237         for(re ll i = 1; i <= N; ++i)
    238         {
    239             std::cin >> xy[i].x;
    240             std::cin >> xy[i].y;
    241             //printf("(%f,%f)
    ", xy[i].x, xy[i].y);
    242         }
    243         printf("Problem %lld:
    ", k);
    244         for(re ll i = 1; i <= M; ++i)
    245         {
    246             ld x, y;
    247             std::cin >> x;
    248             std::cin >> y;
    249             //printf("(%f,%f)
    ", x, y);
    250             if(inpolygon(x, y, N))
    251             {
    252                 printf("Within
    ");
    253             }
    254             else
    255             {
    256                 printf("Outside
    ");
    257             }
    258         }
    259     }
    260     return 0;
    261     */
    262 }

    2. 圆

      1 /*****************************************************************
      2   3 【注意】三角与反三角最多用一次(损失精度大)。
      4 *****************************************************************/
      5 
      6 #include <bits/stdc++.h>
      7 #include <iomanip>
      8 #include <stack>
      9 #include <cstdio>
     10 #include <cmath>
     11 #include <algorithm>
     12 #define re register
     13 #define il inline
     14 #define ll long long
     15 #define ld long double
     16 using namespace std;
     17 const ld PI = 3.14159265358979323846;
     18 const ll MAXN = 5e5+5;
     19 const ld INF = 1e9;
     20 const ld EPS = 1e-10;
     21 
     22 //点坐标
     23 struct POINT
     24 {
     25     ld x, y;
     26     POINT() : x(0), y(0) {}
     27     POINT(ld _x, ld _y) : x(_x), y(_y) {}
     28 };
     29 //向量
     30 typedef POINT VECTOR;
     31 
     32 //符号函数
     33 ll sgn(ld x)    {return x < -EPS ? -1 : x > EPS;}
     34 //向量+向量=向量,点+向量=点
     35 VECTOR operator + (VECTOR a, VECTOR b)  {return {a.x+b.x,a.y+b.y};}
     36 //向量-向量=向量,点-点=向量
     37 VECTOR operator - (VECTOR a, VECTOR b)  {return {a.x-b.x,a.y-b.y};}
     38 //向量*数=向量
     39 VECTOR operator * (VECTOR a, ld k)  {return {a.x*k,a.y*k};}
     40 VECTOR operator * (ld k, VECTOR a)  {return {a.x*k,a.y*k};}
     41 //向量/数=向量
     42 VECTOR operator / (VECTOR a, ld k)  {return {a.x/k,a.y/k};}
     43 //向量==向量,点==点
     44 bool operator == (const VECTOR a, const VECTOR b)   {return !sgn(a.x-b.x) && !sgn(a.y-b.y);}
     45 //向量偏序关系(先x轴再y轴)
     46 bool operator < (const VECTOR a, const VECTOR b)    {return !sgn(a.x-b.x) ? a.y < b.y : a.x < b.x;}
     47 //取下端点
     48 POINT min(POINT p1, POINT p2)   {return p1.y<p2.y ? p1:p2;}
     49 //取上端点
     50 POINT max(POINT p1, POINT p2)   {return p1.y<p2.y ? p2:p1;}
     51 //点积
     52 ld dot(POINT p1, POINT p2, POINT p) {return (p.x-p1.x)*(p2.x-p1.x)+(p.y-p1.y)*(p2.y-p1.y);} //(三点式:共p1端点)
     53 ld dot(VECTOR a, VECTOR b)      {return a.x*b.x+a.y*b.y;}   //向量式
     54 //叉积
     55 ld cross(POINT p1, POINT p2, POINT p)   {return (p.x-p1.x)*(p2.y-p1.y)-(p.y-p1.y)*(p2.x-p1.x);} //(三点式:共p1端点)
     56 ld cross(VECTOR a, VECTOR b)    {return a.x*b.y-a.y*b.x;}   //向量式(aXb)
     57 //向量长度
     58 ld vlen(VECTOR a) {return sqrt(dot(a,a));}
     59 //向量旋转(逆时针)
     60 VECTOR rot(VECTOR a, ld sita)   {return {a.x*cos(sita)-a.y*sin(sita),a.x*sin(sita)+a.y*cos(sita)};}
     61 //取单位向量
     62 VECTOR unit(VECTOR a)   {return a/vlen(a);}
     63 //取正交单位向量
     64 VECTOR norm(VECTOR a)   {return VECTOR(-a.y,a.x)/vlen(a);}
     65 //向量夹角余弦值
     66 ld vcos(VECTOR a, VECTOR b) {return dot(a,b)/(vlen(a)*vlen(b));}
     67 //向量夹角正弦值
     68 ld vsin(VECTOR a, VECTOR b) {return fabs(cross(a,b))/(vlen(a)*vlen(b));}
     69 //向量极角
     70 //【注意】使用系统反三角函数返回nan是因为超出了定义域(大多数由于精度造成)
     71 //【注意】故使用系统反三角函数需处理边界
     72 ld vagl(VECTOR a, VECTOR b)
     73 {
     74     ld d = dot(a,b)/(vlen(a)*vlen(b));
     75     if(sgn(d-1) == 0)   return 0;
     76     else if(sgn(d+1) == 0)  return PI;
     77     else    return acos(d);
     78 }  //向量间[0,PI]
     79 ld vagl(VECTOR a) {return atan2(a.y, a.x);}                 //单向量[-PI,PI]
     80 
     81 //直线
     82 struct LINE
     83 {
     84     POINT p;   //直线上定点
     85     VECTOR v;   //直线方向向量
     86     LINE() {}
     87     //点向式
     88     LINE(POINT _p, VECTOR _v):p(_p), v(_v) {}
     89     //点角式
     90     LINE(POINT _p, ld a):p(_p), v(VECTOR(cos(a), sin(a))) {}
     91     //通过参数t唯一确定直线上点坐标(点向式:P = p + t*v)
     92     POINT point(ld t)   {return p + t*v;}
     93     //平移向量r
     94     LINE trans(VECTOR r)    {return LINE(p+r,v);}
     95 };
     96 
     97 //点到直线距离
     98 ld PLdist(POINT p, LINE L) {return fabs(sin(vagl(p-L.p,L.v))*vlen(p-L.p));}
     99 
    100 //
    101 struct CIRCLE
    102 {
    103     POINT c;        //圆心坐标
    104     ld r;           //半径
    105     CIRCLE() {}
    106     CIRCLE(POINT _c, ld _r):c(_c),r(_r) {}
    107     //通过圆心角唯一确定圆上点坐标
    108     POINT point(ld a)   //a为圆心角
    109     {
    110         return POINT(c.x+cos(a)*r, c.y+sin(a)*r);
    111     }
    112 };
    113 
    114 //两点中垂线
    115 //返回中垂线的个数
    116 ll PerpendicularBisector(POINT p1, POINT p2, LINE &L)
    117 {
    118     if(p1 == p2)    return 0;
    119     L.p = (p1+p2)/2;
    120     L.v = rot(p2-p1, PI/2);
    121     return 1;
    122 }
    123 
    124 //两直线交点
    125 POINT LLitesct;
    126 
    127 //直线与直线的交点
    128 //返回交点个数
    129 ll LineLineIntersection(LINE L1, LINE L2, POINT &p)
    130 {
    131     if(sgn(cross(L1.v, L2.v)) == 0)
    132     {
    133         if(sgn(cross(L2.p-L1.p,L2.v) == 0))    return -1;          //两直线重合
    134         else    return 0;                       //两直线平行
    135     }
    136     //以L1为准
    137     ld t1 = -cross(L2.p-L1.p, L2.v)/cross(L2.v, L1.v);
    138     p = L1.p + t1*L1.v;
    139     //以L2为准
    140     //ld t2 = -cross(L1.p-L2.p, L1.v)/cross(L1.v, L2.v);
    141     //p = L2.p + t2*L2.v;
    142     return 1;
    143 }
    144 
    145 //直线与圆交点数组
    146 vector <POINT> LCitesct;
    147 
    148 //直线与圆的交点
    149 //解方程:(at+b)^2+(ct+d)^2 = r^2;
    150 //化简得:et^2+ft+g = 0
    151 //返回交点个数
    152 ll LineCircleIntersection(LINE L, CIRCLE C, ld &t1, ld &t2, vector <POINT>& sol)
    153 {
    154     ld a = L.v.x, b = L.p.x-C.c.x, c = L.v.y, d = L.p.y-C.c.y;
    155     ld e = a*a + c*c, f = 2*(a*b + c*d), g = b*b + d*d - C.r*C.r;
    156     ld delta = f*f - 4*e*g;
    157     if(sgn(delta) < 0)          //相离
    158         return 0;
    159     else if(sgn(delta) == 0)    //相切
    160     {
    161         t1 = t2 = -f/(2*e);
    162         sol.push_back(L.point(t1));
    163         return 1;
    164     }
    165     else
    166     {
    167         t1 = (-f - sqrt(delta))/(2*e);
    168         t2 = (-f + sqrt(delta))/(2*e);
    169         sol.push_back(L.point(t1));
    170         sol.push_back(L.point(t2));
    171         return 2;
    172     }
    173 }
    174 ll LineCircleIntersection(LINE L, CIRCLE C, POINT *p)
    175 {
    176     ld t1, t2;
    177     ld a = L.v.x, b = L.p.x-C.c.x, c = L.v.y, d = L.p.y-C.c.y;
    178     ld e = a*a + c*c, f = 2*(a*b + c*d), g = b*b + d*d - C.r*C.r;
    179     ld delta = f*f - 4*e*g;
    180     if(sgn(delta) < 0)          //相离
    181         return 0;
    182     else if(sgn(delta) == 0)    //相切
    183     {
    184         t1 = t2 = -f/(2*e);
    185         p[0] = L.point(t1);
    186         return 1;
    187     }
    188     else
    189     {
    190         t1 = (-f - sqrt(delta))/(2*e);
    191         t2 = (-f + sqrt(delta))/(2*e);
    192         p[0] = L.point(t1);
    193         p[1] = L.point(t2);
    194         return 2;
    195     }
    196 }
    197 
    198 //圆与圆的交点数组
    199 vector <POINT> CCitesct;
    200 
    201 //两圆交点
    202 ll CircleCircleIntersection(CIRCLE C1, CIRCLE C2, vector <POINT>& sol)
    203 {
    204     ld d = vlen(C1.c-C2.c);
    205     if(sgn(d) == 0)
    206     {
    207         if(sgn(C1.r-C2.r) == 0) return -1;          //两圆重合
    208         else    return 0;                           //同心相含
    209     }
    210     if(sgn(C1.r+C2.r-d) < 0)    return 0;           //相离
    211     if(sgn(fabs(C1.r-C2.r)-d) > 0)  return 0;       //异心相含
    212 
    213     //设两圆交点为p1,p2
    214     ld a = vagl(C2.c - C1.c);                       //向量C1C2的极角
    215     ld b = acos((C1.r*C1.r + d*d - C2.r*C2.r)/(2*C1.r*d));  //向量C1C2到C1P1(或C1P2)的夹角
    216     POINT p1 = C1.point(a-b);
    217     POINT p2 = C1.point(a+b);
    218     sol.push_back(p1);
    219     if(p1 == p2)    return 1;                       //相切
    220     sol.push_back(p2);                              //相交
    221     return 2;
    222 }
    223 
    224 //圆切线方向向量
    225 VECTOR LCtgt[2];
    226 
    227 //过定点作圆的切线
    228 ll Tangents(POINT p, CIRCLE C, VECTOR *v)
    229 {
    230     VECTOR u = C.c - p;
    231     ld dist = vlen(u);
    232     if(dist < C.r)  return 0;       //点在圆内,没有切线
    233     else if(sgn(dist-C.r) == 0)     //点在圆上,一条切线
    234     {
    235         v[0] = rot(u, PI/2);
    236         return 1;
    237     }
    238     else                            //点在圆外,两条切线
    239     {
    240         ld a = asin(C.r/dist);
    241         v[0] = rot(u, -a);
    242         v[1] = rot(u, +a);
    243         return 2;
    244     }
    245 }
    246 
    247 //两圆公切点(分别在两圆上)
    248 POINT CCtgt1[4], CCtgt2[4];
    249 
    250 //两圆公切线
    251 ll Tangents(CIRCLE C1, CIRCLE C2, POINT *a, POINT *b)
    252 {
    253     ll cnt = 0;
    254     if(C1.r < C2.r) {swap(C1, C2); swap(a, b);}
    255     ld d = dot(C1.c-C2.c, C1.c-C2.c);
    256     ld rdiff = C1.r - C2.r;
    257     ld rsum = C1.r + C2.r;
    258     if(d < rdiff*rsum)  return 0;                       //相含
    259     ld base = atan2(C2.c.y-C1.c.y, C2.c.x-C1.c.x);      //向量C1C2极角[-PI,PI]
    260     if(sgn(d) == 0 && sgn(C1.r-C2.r) == 0)  return -1;  //重合
    261     if(sgn(d-rdiff*rdiff) == 0)                         //内切
    262     {
    263         a[cnt] = C1.point(base);
    264         b[cnt++] = C2.point(base);
    265         return 1;
    266     }
    267     ld ang = acos((C1.r-C2.r)/sqrt(d));                 //夹角[0,PI]
    268     a[cnt] = C1.point(base+ang);
    269     b[cnt++] = C2.point(base+ang);
    270     a[cnt] = C1.point(base-ang);
    271     b[cnt++] = C2.point(base-ang);
    272     if(sgn(d-rsum*rsum))                                //外切
    273     {
    274         a[cnt] = C1.point(base);
    275         b[cnt++] = C2.point(base+PI);
    276     }
    277     else if(d > rsum*rsum)                              //相离
    278     {
    279         ld ang = acos((C1.r+C2.r)/sqrt(d));
    280         a[cnt] = C1.point(base+ang);
    281         b[cnt++] = C2.point(base+ang);
    282         a[cnt] = C1.point(base-ang);
    283         b[cnt++] = C2.point(base-ang);
    284     }
    285     return cnt;                                         //切点个数
    286 }
    287 
    288 //三点外接圆
    289 //返回外接圆的个数
    290 ll CircumscribedCircle(POINT a, POINT b, POINT c, CIRCLE &C)
    291 {
    292     if(sgn(cross(a-b, b-c)) == 0)   return 0;           //三点共线
    293     LINE L1, L2;
    294     PerpendicularBisector(a, b, L1);
    295     PerpendicularBisector(b, c, L2);
    296     POINT p;
    297     LineLineIntersection(L1, L2, p);
    298     C.c = p;
    299     C.r = vlen(p-a);
    300     return 1;
    301 }
    302 
    303 //三点内接圆
    304 //返回内接圆的个数
    305 ll InscribedCircle(POINT a, POINT b, POINT c, CIRCLE &C)
    306 {
    307     if(sgn(cross(a-b, b-c)) == 0)   return 0;           //三点共线
    308     LINE L1 = LINE(a, (vagl(b-a)+vagl(c-a))/2), L2 = LINE(b, (vagl(a-b)+vagl(c-b))/2);
    309     POINT p;
    310     LineLineIntersection(L1, L2, p);
    311     C.c = p;
    312     C.r = PLdist(p,LINE(a,b-a));
    313     return 1;
    314 }
    315 
    316 
    317 int main()
    318 {
    319     //UVA-12304
    320     //ios::sync_with_stdio(false);
    321     string str;
    322     while(cin >> str)
    323     {
    324         if(str == "CircumscribedCircle")
    325         {
    326             POINT a, b, c;
    327             cin >> a.x >> a.y >> b.x >> b.y >> c.x >> c.y;
    328             CIRCLE C;
    329             CircumscribedCircle(a,b,c,C);
    330             cout << fixed << "(" << C.c.x << "," << C.c.y << "," << C.r << ")" << endl;
    331         }
    332         else if(str == "InscribedCircle")
    333         {
    334             POINT a, b, c;
    335             cin >> a.x >> a.y >> b.x >> b.y >> c.x >> c.y;
    336             CIRCLE C;
    337             InscribedCircle(a,b,c,C);
    338             cout << fixed << "(" << C.c.x << "," << C.c.y << "," << C.r << ")" << endl;
    339         }
    340         //求已知圆的切线的极角
    341         else if(str == "TangentLineThroughPoint")
    342         {
    343             CIRCLE C;
    344             POINT p;
    345             cin >> C.c.x >> C.c.y >> C.r >> p.x >> p.y;
    346             VECTOR v[2];
    347             ll cnt = Tangents(p,C,v);
    348             cout << "[";
    349             if(cnt == 1)
    350             {
    351                 ld r = vagl(v[0]);
    352                 ld a = (sgn(r) < 0 ? PI+r : sgn(r-PI) == 0 ? 0 : r)*180/PI;
    353                 cout << fixed << a;
    354             }
    355             else if(cnt == 2)
    356             {
    357                 ld mxr = vagl(v[0]);
    358                 ld mir = vagl(v[1]);
    359                 ld mxa = (sgn(mxr) < 0 ? PI+mxr : sgn(mxr-PI) == 0 ? 0 : mxr)*180/PI;
    360                 ld mia = (sgn(mir) < 0 ? PI+mir : sgn(mir-PI) == 0 ? 0 : mir)*180/PI;
    361                 if(mxa < mia)   swap(mxa,mia);
    362                 cout << fixed << mia << "," << mxa;
    363             }
    364             cout << "]" << endl;
    365         }
    366         //求过定点且与已知直线相切的圆的圆心
    367         else if(str == "CircleThroughAPointAndTangentToALineWithRadius")
    368         {
    369             CIRCLE P;
    370             POINT p1, p2;
    371             cin >> P.c.x >> P.c.y >> p1.x >> p1.y >> p2.x >> p2.y >> P.r;
    372             VECTOR v = rot(p2-p1,PI/2);
    373             ld t = P.r/vlen(v);
    374             LINE L1 = LINE(p1+t*v,p2-p1), L2 = LINE(p1-t*v,p2-p1);
    375             ll cnt = 0;
    376             ld t1, t2;
    377             vector <POINT> C;
    378             cnt += LineCircleIntersection(L1,P,t1,t2,C);
    379             cnt += LineCircleIntersection(L2,P,t1,t2,C);
    380             sort(C.begin(),C.end());
    381             cout << "[";
    382             if(cnt == 1)    cout << fixed << "(" << C[0].x << "," << C[0].y << ")";
    383             else if(cnt == 2)   cout << fixed << "(" << C[0].x << "," << C[0].y << ")" << "," << "(" << C[1].x << "," << C[1].y << ")";
    384             cout << "]" << endl;
    385         }
    386         //求给定半径且与两条不相交直线相切的圆的圆心
    387         else if(str == "CircleTangentToTwoLinesWithRadius")
    388         {
    389             POINT p1, p2, p3, p4;
    390             ld r;
    391             cin >> p1.x >> p1.y >> p2.x >> p2.y >> p3.x >> p3.y >> p4.x >> p4.y >> r;
    392             LINE L1 = LINE(p1,p2-p1), L2 = LINE(p3,p4-p3);
    393             VECTOR v1 = norm(p2-p1), v2 = norm(p4-p3);
    394             LINE L11 = L1.trans(r*v1), L12 = L1.trans(-r*v1);
    395             LINE L21 = L2.trans(r*v2), L22 = L2.trans(-r*v2);
    396             vector <POINT> C;
    397             POINT p;
    398             LineLineIntersection(L11,L21,p);
    399             C.push_back(p);
    400             LineLineIntersection(L11,L22,p);
    401             C.push_back(p);
    402             LineLineIntersection(L12,L21,p);
    403             C.push_back(p);
    404             LineLineIntersection(L12,L22,p);
    405             C.push_back(p);
    406             sort(C.begin(),C.end());
    407             cout << "[";
    408             cout << fixed << "(" << C[0].x << "," << C[0].y << ")" << "," << "(" << C[1].x << "," << C[1].y << ")" << ",";
    409             cout << fixed << "(" << C[2].x << "," << C[2].y << ")" << "," << "(" << C[3].x << "," << C[3].y << ")";
    410             cout << "]" << endl;
    411 
    412         }
    413         else if(str == "CircleTangentToTwoDisjointCirclesWithRadius")
    414         {
    415             CIRCLE C1, C2;
    416             ld r;
    417             cin >> C1.c.x >> C1.c.y >> C1.r >> C2.c.x >> C2.c.y >> C2.r >> r;
    418             C1.r += r, C2.r += r;
    419             vector <POINT> p;
    420             ll cnt = CircleCircleIntersection(C1,C2,p);
    421             sort(p.begin(),p.end());
    422             cout << "[";
    423             if(cnt == 1)
    424             {
    425                 cout << fixed << "(" << p[0].x << "," << p[0].y << ")";
    426             }
    427             else if(cnt == 2)
    428             {
    429                 cout << fixed << "(" << p[0].x << "," << p[0].y << ")" << "," << "(" << p[1].x << "," << p[1].y << ")";
    430             }
    431             cout << "]" << endl;
    432         }
    433     }
    434     return 0;
    435 }

    3. 半平面

      1 /*****************************************************************
      2                             半平面交
      3 【注意】三角与反三角最多用一次(损失精度大)。
      4 【注意】数组下标一律从1开始。
      5 *****************************************************************/
      6 
      7 //#include <bits/stdc++.h>
      8 #include <iostream>
      9 #include <cstring>
     10 #include <cstdlib>
     11 #include <iomanip>
     12 #include <stack>
     13 #include <cstdio>
     14 #include <cmath>
     15 #include <algorithm>
     16 #define re register
     17 #define il inline
     18 #define ll int
     19 #define ld double
     20 using namespace std;
     21 const ld PI = 3.14159265358979323846;
     22 const ll MAXN = 5e5+5;
     23 const ld INF = 1e9;
     24 const ld EPS = 1e-8;
     25 
     26 //点坐标
     27 struct POINT
     28 {
     29     ld x, y;
     30     POINT() : x(0), y(0) {}
     31     POINT(ld _x, ld _y) : x(_x), y(_y) {}
     32 };
     33 //向量
     34 typedef POINT VECTOR;
     35 
     36 //符号函数
     37 ll sgn(ld x)    {return x < -EPS ? -1 : x > EPS;}
     38 //向量+向量=向量,点+向量=点
     39 VECTOR operator + (VECTOR a, VECTOR b)  {return {a.x+b.x,a.y+b.y};}
     40 //向量-向量=向量,点-点=向量
     41 VECTOR operator - (VECTOR a, VECTOR b)  {return {a.x-b.x,a.y-b.y};}
     42 //向量*数=向量
     43 VECTOR operator * (VECTOR a, ld k)  {return {a.x*k,a.y*k};}
     44 VECTOR operator * (ld k, VECTOR a)  {return {a.x*k,a.y*k};}
     45 //向量/数=向量
     46 VECTOR operator / (VECTOR a, ld k)  {return {a.x/k,a.y/k};}
     47 //向量==向量,点==点
     48 bool operator == (const VECTOR a, const VECTOR b)   {return !sgn(a.x-b.x) && !sgn(a.y-b.y);}
     49 //向量偏序关系(先x轴再y轴)
     50 bool operator < (const VECTOR a, const VECTOR b)    {return !sgn(a.x-b.x) ? a.y < b.y : a.x < b.x;}
     51 //取下端点
     52 POINT min(POINT p1, POINT p2)   {return p1.y<p2.y ? p1:p2;}
     53 //取上端点
     54 POINT max(POINT p1, POINT p2)   {return p1.y<p2.y ? p2:p1;}
     55 //点积
     56 ld dot(POINT p1, POINT p2, POINT p) {return (p.x-p1.x)*(p2.x-p1.x)+(p.y-p1.y)*(p2.y-p1.y);} //(三点式:共p1端点)
     57 ld dot(VECTOR a, VECTOR b)      {return a.x*b.x+a.y*b.y;}   //向量式
     58 //叉积
     59 ld cross(POINT p1, POINT p2, POINT p)   {return (p.x-p1.x)*(p2.y-p1.y)-(p.y-p1.y)*(p2.x-p1.x);} //(三点式:共p1端点)
     60 ld cross(VECTOR a, VECTOR b)    {return a.x*b.y-a.y*b.x;}   //向量式(aXb)
     61 //向量长度
     62 ld vlen(VECTOR a) {return sqrt(dot(a,a));}
     63 //向量旋转(逆时针)
     64 VECTOR rot(VECTOR a, ld sita)   {return {a.x*cos(sita)-a.y*sin(sita),a.x*sin(sita)+a.y*cos(sita)};}
     65 //取单位向量
     66 VECTOR unit(VECTOR a)   {return a/vlen(a);}
     67 //取正交单位向量
     68 VECTOR norm(VECTOR a)   {return VECTOR(-a.y,a.x)/vlen(a);}
     69 //向量夹角余弦值
     70 ld vcos(VECTOR a, VECTOR b) {return dot(a,b)/(vlen(a)*vlen(b));}
     71 //向量夹角正弦值
     72 ld vsin(VECTOR a, VECTOR b) {return fabs(cross(a,b))/(vlen(a)*vlen(b));}
     73 //向量极角
     74 //【注意】使用系统反三角函数返回nan是因为超出了定义域(大多数由于精度造成)
     75 //【注意】故使用系统反三角函数需处理边界
     76 ld vagl(VECTOR a, VECTOR b)
     77 {
     78     ld d = dot(a,b)/(vlen(a)*vlen(b));
     79     if(sgn(d-1) == 0)   return 0;
     80     else if(sgn(d+1) == 0)  return PI;
     81     else    return acos(d);
     82 }  //向量间[0,PI]
     83 ld vagl(VECTOR a) {return atan2(a.y, a.x);}                 //单向量[-PI,PI]
     84 
     85 //凸包
     86 ll vexcnt = 0;      //凸包顶点数
     87 POINT xy[MAXN], xyt[MAXN];     //多边形顶点存储,临时多边形顶点存储
     88 POINT convex[MAXN]; //凸包顶点数组
     89 
     90 //求凸包周长
     91 ld convexperimeter(POINT *poly, ll n)
     92 {
     93     ld ans = 0;
     94     for(re ll i = 1; i <= n; ++i)
     95     {
     96         ans += vlen(poly[i]-poly[i%n+1]);
     97     }
     98     return ans;
     99 }
    100 //计算多边形面积
    101 ld polygonarea(POINT *poly, ll n)
    102 {
    103     ld aera = 0;
    104     for(re ll i = 0; i < n; ++i)
    105     {
    106         aera += cross(poly[i+1], poly[(i+1)%n+1]);
    107     }
    108     return fabs(aera)/2;        //不加绝对值为有向面积
    109 }
    110 
    111 //andrew算法求凸包(凸包数组正序为逆时针)
    112 //1.严格凸包(没有重复点和三点共线)2.非严格凸包(允许重复点和三点共线)
    113 void andrew(ll n)
    114 {
    115     vexcnt = 0;                            //初始化数据
    116     memcpy(xyt,xy,sizeof(xy));          //备份原来顶点集
    117     sort(xyt+1,xyt+n+1);                //排序后xyt[1]和xyt[n]一定在凸包中
    118     //求解凸包顶点集
    119     //正向扫描(下凸包)
    120     convex[++vexcnt] = xyt[1];          //xyt[1]一定在凸包中
    121     ll j = 2;
    122     while(j <= n)
    123     {
    124         if(vexcnt <= 1)                 //因为xyt[2]不一定在凸包集中
    125         {
    126             convex[++vexcnt] = xyt[j++];
    127         }
    128         else
    129         {
    130             //取前面两个点
    131             //严格凸包:sgn(cross(convex[vexcnt]-convex[vexcnt-1],xyt[j]-convex[vexcnt]))>0
    132             //非严格凸包:sgn(cross(convex[vexcnt]-convex[vexcnt-1],xyt[j]-convex[vexcnt]))>=0
    133             if(sgn(cross(convex[vexcnt]-convex[vexcnt-1],xyt[j]-convex[vexcnt]))>=0)  convex[++vexcnt] = xyt[j++];
    134             else    --vexcnt;
    135         }
    136     }
    137     //反向扫描(上凸包)
    138     //至少包含了xyt[1]和xyt[n]
    139     ll record = vexcnt;                 //记录下凸包的结果
    140     ll k = n-1;
    141     while(k >= 1)
    142     {
    143         if(vexcnt <= record)
    144         {
    145             convex[++vexcnt] = xyt[k--];
    146         }
    147         else
    148         {
    149             //取前面两个点
    150             //严格凸包:sgn(cross(convex[vexcnt]-convex[vexcnt-1],xyt[k]-convex[vexcnt]))>0
    151             //非严格凸包:sgn(cross(convex[vexcnt]-convex[vexcnt-1],xyt[k]-convex[vexcnt]))>=0
    152             if(sgn(cross(convex[vexcnt]-convex[vexcnt-1],xyt[k]-convex[vexcnt]))>0)  convex[++vexcnt] = xyt[k--];
    153             else    --vexcnt;
    154         }
    155     }
    156     while(convex[--vexcnt]==convex[1]);     //因为和xyt[1]相等的点都在下凸包中了
    157     //for(re ll i = 1; i <= vexcnt; ++i)    cout << "(" << convex[i].x << "," << convex[i].y << ")" << endl;
    158     return;
    159 }
    160 
    161 //直线
    162 struct LINE
    163 {
    164     POINT p;   //直线上定点
    165     VECTOR v;   //直线方向向量
    166     LINE() {}
    167     //点向式
    168     LINE(POINT _p, VECTOR _v):p(_p), v(_v) {}
    169     //点角式
    170     LINE(POINT _p, ld a):p(_p), v(VECTOR(cos(a), sin(a))) {}
    171     //通过参数t唯一确定直线上点坐标(点向式:P = p + t*v)
    172     POINT point(ld t)   {return p + t*v;}
    173     //平移向量r
    174     LINE trans(VECTOR r)    {return LINE(p+r,v);}
    175 };
    176 
    177 //有向直线
    178 struct DIRLINE
    179 {
    180     POINT p;        //定点
    181     VECTOR v;       //方向向量
    182     ld a;           //极角
    183     DIRLINE() {}
    184     DIRLINE(POINT _p, VECTOR _v):p(_p),v(_v),a(atan2(v.y,v.x)) {}
    185     bool operator < (const DIRLINE &DL) const    {return a < DL.a;}
    186 };
    187 
    188 //点到直线距离
    189 ld PLdist(POINT p, LINE L) {return fabs(sin(vagl(p-L.p,L.v))*vlen(p-L.p));}
    190 //点在有向直线左边
    191 bool OnLeft(DIRLINE DL, POINT p) {return sgn(cross(DL.v,p-DL.p)) > 0;}
    192 //点在有向直线右边
    193 bool OnRight(DIRLINE DL, POINT p)   {return sgn(cross(DL.v,p-DL.p)) < 0;}
    194 
    195 //两点中垂线
    196 //返回中垂线的个数
    197 ll PerpendicularBisector(POINT p1, POINT p2, LINE &L)
    198 {
    199     if(p1 == p2)    return 0;
    200     L.p = (p1+p2)/2;
    201     L.v = rot(p2-p1, PI/2);
    202     return 1;
    203 }
    204 
    205 //两直线交点
    206 POINT LLitesct;
    207 
    208 //直线与直线的交点
    209 //返回交点个数
    210 ll LineLineIntersection(LINE L1, LINE L2, POINT &p)
    211 {
    212     if(sgn(cross(L1.v, L2.v)) == 0)
    213     {
    214         if(sgn(cross(L2.p-L1.p,L2.v)) == 0)   return -1;          //两直线重合
    215         else    return 0;                       //两直线平行
    216     }
    217     //以L1为准
    218     ld t1 = -cross(L2.p-L1.p, L2.v)/cross(L2.v, L1.v);
    219     p = L1.p + t1*L1.v;
    220     //以L2为准
    221     //ld t2 = -cross(L1.p-L2.p, L1.v)/cross(L1.v, L2.v);
    222     //p = L2.p + t2*L2.v;
    223     return 1;
    224 }
    225 ll LineLineIntersection(DIRLINE DL1, DIRLINE DL2, POINT &p)
    226 {
    227     if(sgn(cross(DL1.v, DL2.v)) == 0)
    228     {
    229         if(sgn(cross(DL2.p-DL1.p,DL2.v)) == 0)    return -1;          //两直线重合
    230         else    return 0;                       //两直线平行
    231     }
    232     //以L1为准
    233     ld t1 = -cross(DL2.p-DL1.p, DL2.v)/cross(DL2.v, DL1.v);
    234     p = DL1.p + t1*DL1.v;
    235     //以L2为准
    236     //ld t2 = -cross(L1.p-L2.p, L1.v)/cross(L1.v, L2.v);
    237     //p = L2.p + t2*L2.v;
    238     return 1;
    239 }
    240 
    241 //半平面交
    242 //返回凸包顶点数
    243 //【注意】若半平面交可能是无界区域,则需在外面加一个坐标很大的包围框
    244 ll HalfPlaneIntersection(DIRLINE *DL, ll n, POINT *poly)
    245 {
    246     sort(DL+1,DL+n+1);
    247     ll first, last;                 //双端队列的首尾元素下标
    248     DIRLINE *q = new DIRLINE[n+1];    //双端队列
    249     POINT *p = new POINT[n+1];        //p[i]为q[i]和q[i+1]的交点
    250     first = last = 1;
    251     q[1] = DL[1];                   //初始化只有一个半平面L[0]
    252     for(re ll i = 2; i <= n; ++i)
    253     {
    254         while(first < last && !OnLeft(DL[i], p[last-1])) --last;
    255         while(first < last && !OnLeft(DL[i], p[first])) ++first;
    256         q[++last] = DL[i];
    257         if(sgn(cross(q[last].v, q[last-1].v)) == 0)         //两向量平行
    258         {
    259             --last;
    260             if(OnLeft(q[last], DL[i].p)) q[last] = DL[i];   //同向取内测的一个
    261         }
    262         if(first < last)    LineLineIntersection(q[last-1],q[last],p[last-1]);
    263     }
    264     while(first < last && !OnLeft(q[first],p[last-1]))  --last;
    265     //删除无用平面
    266     if(last-first <= 1) return 0;                       //空集
    267     LineLineIntersection(q[last],q[first],p[last]);     //计算首尾两个半平面交的交点
    268     ll cnt = 0;
    269     for(re ll i = first; i <= last; ++i)    poly[++cnt] = p[i];
    270     return cnt;
    271 }
    272 
    273 int main()
    274 {
    275     //UVALive-3890
    276     ios::sync_with_stdio(false);
    277     ll n;
    278     //scanf("%lld", &n);
    279     while(cin >> n)
    280     //while(~scanf("%d", &n))
    281     {
    282         if(n == 0)  break;
    283         POINT *p = new POINT[n+1];
    284         for(re ll i = 1; i <= n; ++i)
    285         {
    286             cin >> p[i].x >> p[i].y;
    287             //scanf("%lf%lf%lf%lf", &p1.x, &p1.y, &p2.x, &p2.y);
    288         }
    289         VECTOR *v = new VECTOR[n+1];
    290         VECTOR *r = new VECTOR[n+1];
    291         for(re ll i = 1; i <= n; ++i)
    292         {
    293             v[i] = p[i%n+1]-p[i];
    294             r[i] = norm(v[i]);      //单位正交向量
    295         }
    296         //二分答案
    297         ld left = 0, right = 20000;
    298         DIRLINE *DL = new DIRLINE[n+1];
    299         POINT *poly = new POINT[n+1];
    300         while(right-left > 1e-6)
    301         {
    302             ld mid = (left+right)/2;
    303             for(re ll i = 1; i <= n; ++i)
    304             {
    305                 DL[i] = DIRLINE(p[i]+r[i]*mid,v[i]);
    306             }
    307             ll cnt = HalfPlaneIntersection(DL,n,poly);
    308             if(!cnt)    right = mid;
    309             else    left = mid;
    310         }
    311         cout << fixed << setprecision(6) << left << endl;
    312     }
    313     return 0;
    314 }
  • 相关阅读:
    简便的将DataSet导入到数据库中
    数据类型的小小研究:Access与SQL Server的数据类型
    【jxust acm 20120708】
    【D ECJTU_ACM 11级队员2012年暑假训练赛(2)】
    【hdu 2101 A + B Problem Too】
    【hdu 1014】
    【hdu 1164 Eddy's research I】
    【开始,安全编程】
    【hdu 1285 确定比赛名次】
    【hdu 1163 Eddy's digital Roots 】
  • 原文地址:https://www.cnblogs.com/BlueHeart0621/p/12238221.html
Copyright © 2011-2022 走看看