哈工大研究生课程- 高级算法设计与分析-课程设计.zip

# 实验部分 ## 实验一: 分治算法 ### 实验目的 * 掌握分治算法的设计思想与方法 * 熟练使用高级编程语言实现分治算法 * 通过对比简单算法以及不同的分治求解思想，理解算法复杂度 ### 实验问题 求解凸包问题：输入是平面上 $n$ 个点的集合$Q$，凸包问题是要输出一个$Q$的凸包。其中，$Q$ 的凸包是一个凸多边形$P$，$Q$ 中的点或者在$P$上或者在$P$中。（详情请见课件） ### 实验步骤 #### 实现基于枚举方法的凸包求解算法 * 算法思想  算法需要四层循环，每层遍历所有点，复杂度为$O(n^4)$ * 实现如下，详见文件[BruteForceCH.hpp](lab01/BruteForceCH.hpp) ```c++ /** * time complexity: O(n^4+nlog(n)) */ vector<Point<int>> BruteForceCH::work() { int n = P.size(); if (n == 0) return vector<Point<int>>(); for (int a = 0; a < n; ++a) { if (P[a].x < 0) continue; // Point::x < 0 -> point deleted for (int b = 0; b < n; ++b) { if (P[a].x < 0) break; // important if (b == a || P[b].x < 0) continue; for (int c = 0; c < n; ++c) { if (P[b].x < 0 || P[a].x < 0) break; // important if (c == a || c == b || P[c].x < 0) continue; // check if abc on same line, if so, delete the middle point if (P[a].cross(P[b], P[c]) == 0) { int i = a, j = b, k = c; if (P[j] < P[i]) swap(i, j); if (P[j] < P[k]) P[j].x = -1; else if (P[i] < P[k]) P[k].x = -1; else P[i].x = -1; continue; } for (int d = 0; d < n; ++d) { if (d == a || d == b || d == c || P[d].x < 0) continue; if (P[a].cross(P[b], P[c])*P[a].cross(P[b], P[d]) >= 0 && P[a].cross(P[c], P[b])*P[a].cross(P[c], P[d]) >= 0 && P[b].cross(P[c], P[a])*P[b].cross(P[c], P[d]) >= 0 ) { P[d].x = -1; // delete } } } } } // A, B, left most & right most int m = 1, i = 0; while (P[i].x < 0) ++i; int a = i, b = i; for (++i; i < n; ++i) { if (P[i].x < 0) continue; ++m; if (P[i] < P[a]) a = i; if (P[b] < P[i]) b = i; } vector<Point<int>> CH(m), S_U(m); CH[0] = P[a]; // A int il = 1, iu = 0; // S_L, S_U for (int i = 0; i < n; ++i) { if (P[i].x < 0 || i == a || i == b) continue; if (P[a].cross(P[b], P[i]) < 0) { // < 0, S_L CH[il++] = P[i]; } else { // > 0, S_U S_U[iu++] = P[i]; } } // A, S_L, B, reverse(S_U) sort(CH.begin() + 1, CH.begin() + il); // S_L sort(S_U.begin(), S_U.begin() + iu); CH[il] = P[b]; // B for (int i = 0; i < iu; ++i) CH[il+1+i] = S_U[iu-1-i]; // reverse(S_U) return CH; } ``` * 首先四层循环，不断删除点，如果这个点在某个三角形内的话，删除的点的$x$坐标标记为-1。 * 需要注意和伪代码略有不同的是（伪代码只描述了主要逻辑，具体实现细节有的地方需要注意），三角形的三点如果共线，那么删除掉三点中中间的点，不进入4层循环，而是返回到删除点的那一层（向外可能多次break） * 最后就是图包点的顺序调整（和伪代码一致），返回 #### 实现基于Graham-Scan的凸包求解算法 * 算法思想  算法需要用到排序，每个点最多访问两次（加入，删除），所以最终复杂度为$O(n\log n)$ * 实现如下，详见文件[GrahamScanCH.hpp](lab01/GrahamScanCH.hpp) ```c++ /** * time complexity: O(nlog(n)) */ vector<Point<int>> GrahamScanCH::work() { int n = P.size(); if (n == 0) return vector<Point<int>>(); swap(P[0], *min_element(P.begin(), P.end())); if (n <= 2) return P; sort(P.begin()+1, P.end(), [&](const auto& a, const auto& b) { int c = P[0].cross(a, b); return c > 0 || (c == 0 && P[0].dist2(a) < P[0].dist2(b)); }); P.push_back(P[0]); // 处理最后一个点更简单 vector<Point<int>> CH; int m = 0; for (int i = 0; i <= n; ++i) { while (m >= 2 && CH[m-2].cross(CH[m-1], P[i]) <= 0) { if (i == n && m == 2) break; // 防止所有点都共线，删除了一端端点 CH.pop_back(), --m; } CH.push_back(P[i]); ++m; } CH.pop_back(); return CH; } ``` * 首先寻找左下角的点，按照它为极点将其余点按照极角排序（极角大小的判定通过叉积计算，比算角度方便） * 为了处理最后一个点方便，额外将极点再加入到末尾 * 循环中注意所有点如果在同一条线上的话会多把一个端点删除掉，所以循环中有一行判断 * 最后弹出额外加入的极点，返回凸包结果 #### 实现基于分治思想的凸包求解算法 * 算法思想 * 实现了两种分治算法，分别见[DivAndConCH.hpp](lab01/DivAndConCH.hpp)和[DivAndConCH2.hpp](lab01/DivAndConCH2.hpp) * DivAndConCH2中的算法采用的是课件中的方式，主要在于merge步骤，前面部分都相同。  * DivAndConCH中的算法采用的[Convex Hull: Divide and Conquer](http://web.ntnu.edu.tw/~algo/ConvexHull.html#6)方式，merge过程通过找上下公切线实现，从左凸包的最右点和右凸包的最左点开始，不断旋转端点到下一点，直到卡死，复杂度也是$O(n)$  * 算法实现 * DivAndConCH，找上下公切线方式merge ```c++ /** * time complexity: O(nlog(n)) */ vector<Point<int>> DivAndConCH::work() { int n = P.size(); if (n == 0) return vector<Point<int>>(); sort(P.begin(), P.end()); return _div_con(0, n); } vector<Point<int>> DivAndConCH::_div_con(int l, int r) { int m = r - l; if (m < 3) return vector<Point<int>>(P.begin()+l, P.begin()+r); vector<Point<int>> lCH = _div_con(l, l + m/2), rCH = _div_con(l + m/2, r); int lsz = lCH.size(), rsz = rCH.size(); // merge int a = 0, b = 0; for (int i = 0; i < lsz; ++i) { if (lCH[a] < lCH[i]) a = i; } for (int i = 0; i < rsz; ++i) { if (rCH[i] < rCH[b]) b = i; } int llb = a, lrb = a, rlb = b, rrb = b; // 下公切线 while (true) { bool changed = false; for (int nlrb = (lrb-1+lsz)%lsz, cs = rCH[rlb].cross(lCH[lrb], lCH[nlrb]); cs > 0 || (cs == 0 && rCH[rlb].dist2(lCH[lrb]) < rCH[rlb].dist2(lCH[nlrb])); nlrb = (lrb-1+lsz)%lsz, cs = rCH[rlb].cross(lCH[lrb], lCH[nlrb])) { lrb = nlrb; changed = true; } for (int nrlb = (rlb+1)%rsz, cs = lCH[lrb].cross(rCH[rlb], rCH[nrlb]); cs < 0 || (cs == 0 && lCH[lrb].dist2(rCH[rlb]) < lCH[lrb].dist2(rCH[nrlb])); nrlb = (rlb+1)%rsz, cs = lCH[lrb].cross(rCH[rlb], rCH[nrlb])) { rlb = nrlb; changed = true; } if (!changed) break; } // 上公切线 while (true) { bool changed = false; for (int nllb = (llb+1)%lsz, cs = rCH[rrb].cross(lCH[llb], lCH[nllb]); cs < 0 || (cs == 0 && rCH[rrb].dist2(lCH[llb]) < rCH[rrb].dist2(lCH[nllb]