天气与日历 切换到窄版

 找回密码
 立即注册
中国膜结构网
十大进口膜材评选 十大国产膜材评选 十大膜结构设计评选 十大膜结构公司评选
查看: 180|回复: 5

GJK算法计算凸多边形之间的距离

[复制链接]
  • TA的每日心情
    开心
    昨天 06:36
  • 签到天数: 15 天

    [LV.4]偶尔看看III

    105

    主题

    11

    回帖

    1308

    积分

    管理员

    积分
    1308
    QQ
    发表于 2024-3-16 09:02:01 | 显示全部楼层 |阅读模式
    1. d = c2 - c1; // 和GJK碰撞检测类似
    2. Simplex.add(support(shape1, shape2, d)); // Simplex 中加入 a 点
    3. Simplex.add(support(shape1, shape2, -d));  // Simplex 中加入 b 点
    4. // 从原点指向 ab 线段上距离原点最近的点的向量, 例如恰好就是答案的话, 则 d.magnitude() 就是答案
    5. d = ClosestPointToOrigin(Simplex.a, Simplex.b);
    6. while (true) {
    7.   // 上面的 ClosestPointToOrigin 得到的 d 方向是从原点到 ab 线段上的最近点的, 所以将d反向, 则指向原点
    8.   d.negate();
    9.   if (d.isZero()) {
    10.     // 则 原点在 Minkowski 和中, 所以发生了碰撞, 但是这当然不是碰撞的唯一情形
    11.     return false;
    12.   }
    13.   // 计算出新的点c, 注意, 因为 d 是朝向坐标原点的
    14.   c = support(shape1, shape2, d);
    15.   // c 在 d 方向上的投影
    16.   double dc = c.dot(d);
    17.   // a 在 d 方向上的投影
    18.   double da = Simplex.a.dot(d);
    19.   // 这个 tolerance 其实对于多边形而言, 取 0 就好了, tolerance 是对于 曲边形才需要设置
    20.   if (|dc - da| < tolerance) {
    21.     distance = d.magnitude();
    22.     return true;
    23.   }
    24.   // 因为我们已经知道了 c 比 a、b更接近原点, 所以只需要保留 a、b中更接近原点的那个点了
    25.   // 然后保留下来的点和c共同构成新的单纯形(即一条线段)
    26.   p1 = ClosestPointToOrigin(Simplex.a, c);
    27.   p2 = ClosestPointToOrigin(Simplex.b, c);
    28.   if (p1.magnitude() < p2.magnitude()) {
    29.     Simplex.b = c;
    30.     d = p1;
    31.   } else {
    32.     Simplex.a = c;
    33.     d = p2;
    34.   }
    35. }
    复制代码

     

     

     

     

    GJK算法计算凸多边形之间的距离
  • TA的每日心情
    开心
    昨天 06:36
  • 签到天数: 15 天

    [LV.4]偶尔看看III

    105

    主题

    11

    回帖

    1308

    积分

    管理员

    积分
    1308
    QQ
     楼主| 发表于 2024-3-16 09:02:14 | 显示全部楼层
    1. d = (11.5, 4.0) - (5.5, 8.5) = (6, -4.5)
    2. Simplex.add(support(shape1, shape2, d)) = (9, 9) - (8, 6) = (1, 3)
    3. Simplex.add(support(shape1, shape2, -d)) = (4, 11) - (13, 1) = (-9, 10)
    4. // 计算新的 d, 直接看 Figure 4 即可,因为 ab 线段上距离原点最近的点就是 (1,3), 所以 d = (1, 3)
    5. d = (1, 3)
    6. d = (-1, -3) // d 反向
    复制代码

     

     

     

     

    GJK算法计算凸多边形之间的距离
  • TA的每日心情
    开心
    昨天 06:36
  • 签到天数: 15 天

    [LV.4]偶尔看看III

    105

    主题

    11

    回帖

    1308

    积分

    管理员

    积分
    1308
    QQ
     楼主| 发表于 2024-3-16 09:02:25 | 显示全部楼层
    1. // 用support 函数计算新的点
    2. c = support(shape1, shape2, d) = (4, 5) - (15, 6) = (-11, -1)
    3. // 14 - (-10) = 24 还不够小,  所以我们不能终止循环, 事实上, 对于凸多边形的情况, 只有为 0 了, 才算够小
    4. dc = 11 + 3 = 14
    5. da = -1 - 9 = -10
    6. // 边 AC [(1, 3) to (-11, -1)] 和 边 BC [(-9, 10) to (-11, -1)]相比更加接近原点, 所以应该蒯掉 b 点, 也就是将 b 替换为 新发现的, 更加接近原点的 c
    7. b = c
    8. // 设置新的 d, 我们实际上是朝着原点在不断进发
    9. d = p1
    复制代码

     

     

     

     

    GJK算法计算凸多边形之间的距离
  • TA的每日心情
    开心
    昨天 06:36
  • 签到天数: 15 天

    [LV.4]偶尔看看III

    105

    主题

    11

    回帖

    1308

    积分

    管理员

    积分
    1308
    QQ
     楼主| 发表于 2024-3-16 09:02:33 | 显示全部楼层
    1. d = (1.07, -1.34) // 即上面的 p1 的相反向量 -p1
    2. // 根据
    3. c = support(shape1, shape2, d) = (9, 9) - (8, 6) = (1, 3)
    4. // 够小了!
    5. dc = -1.07 + 4.02 = 2.95
    6. da = -1.07 + 4.02 = 2.95
    7. // ||d|| = 1.7147886 我们完成了迭代
    8. distance = 1.71
    复制代码

     

     

     

     

    GJK算法计算凸多边形之间的距离
  • TA的每日心情
    开心
    昨天 06:36
  • 签到天数: 15 天

    [LV.4]偶尔看看III

    105

    主题

    11

    回帖

    1308

    积分

    管理员

    积分
    1308
    QQ
     楼主| 发表于 2024-3-16 09:02:46 | 显示全部楼层
    1. 题目概述
    2. 给定两个不相交的凸多边形,求其之间最近距离
    3. 时限
    4. 1000ms 64MB
    5. 输入
    6. 第一行正整数N,M,代表两个凸多边形顶点数,其后N行,每行两个浮点数x,y,描述多边形1的一个点的坐标,其后M
    7. 行,每行两个浮点数x,y,描述多边形2的一个点的坐标,输入到N=M=0为止
    8. 输入保证是按照顺时针或者逆时针给出凸包上的点.
    9. 限制
    10. 3<=N,M<=10000;-10000<=x,y<=10000
    11. 输出
    12. 每行一个浮点数,为所求最近距离,误差在1e-3内均视为正确
    13. 样例输入
    14. 4 4
    15. 0.00000 0.00000
    16. 0.00000 1.00000
    17. 1.00000 1.00000
    18. 1.00000 0.00000
    19. 2.00000 0.00000
    20. 2.00000 1.00000
    21. 3.00000 1.00000
    22. 3.00000 0.00000
    23. 0 0
    24. 样例输出
    25. 1.00000
    复制代码

     

     

     

     

    GJK算法计算凸多边形之间的距离
  • TA的每日心情
    开心
    昨天 06:36
  • 签到天数: 15 天

    [LV.4]偶尔看看III

    105

    主题

    11

    回帖

    1308

    积分

    管理员

    积分
    1308
    QQ
     楼主| 发表于 2024-3-16 09:03:45 | 显示全部楼层
    1. //#include "stdafx.h"
    2. //#define LOCAL
    3. #pragma GCC optimize(2)
    4. #pragma G++ optimize(2)
    5. #pragma warning(disable:4996)
    6. //#pragma comment(linker, "/STACK:1024000000,1024000000")
    7. #include <stdio.h>
    8. #include <iostream>
    9. #include <iomanip>
    10. #include <string>
    11. #include <ctype.h>
    12. #include <string.h>
    13. #include <fstream>
    14. #include <sstream>
    15. #include <math.h>
    16. #include <map>
    17. //#include <unordered采用map>
    18. #include <algorithm>
    19. #include <vector>
    20. #include <stack>
    21. #include <queue>
    22. #include <set>
    23. #include <time.h>
    24. #include <stdlib.h>
    25. #include <bitset>
    26. using namespace std;
    27. //#define int unsigned long long
    28. //#define int long long
    29. #define re register int
    30. #define ci const int
    31. #define ui unsigned int
    32. typedef pair<int, int> P;
    33. #define FE(cur) for(re h = head[cur], to; ~h; h = g[h].nxt)
    34. #define ilv inline void
    35. #define ili inline int
    36. #define ilc inline char
    37. #define ild inline double
    38. #define ilp inline P
    39. #define LEN(cur) (hjt[cur].r - hjt[cur].l)
    40. #define MID(cur) (hjt[cur].l + hjt[cur].r >> 1)
    41. #define SQUARE(x) ((x) * (x))
    42. typedef vector<int>::iterator vit;
    43. typedef set<int>::iterator sit;
    44. typedef map<int, int>::iterator mit;
    45. const int inf = ~0u>>1;
    46. const double PI = acos(-1.0), eps = 1e-8;
    47. namespace fastio
    48. {
    49.     const int BUF = 1 << 21;
    50.     char fr[BUF], fw[BUF], *pr1 = fr, *pr2 = fr;int pw;
    51.     ilc gc() { return pr1 == pr2 && (pr2 = (pr1 = fr) + fread(fr, 1, BUF, stdin), *pr2 = 0, pr1 == pr2) ? EOF : *pr1++; }
    52.     ilv flush() { fwrite(fw, 1, pw, stdout); pw = 0; }
    53.     ilv pc(char c) { if (pw >= BUF) flush(); fw[pw++] = c; }
    54.     ili read(int &x)
    55.     {
    56.         x = 0; int f = 1; char c = gc(); if (!~c) return EOF;
    57.         while(!isdigit(c)) { if (c == '-') f = -1; c = gc(); }
    58.         while(isdigit(c)) x = (x << 3) + (x << 1) + (c ^ 48), c = gc();
    59.         x *= f; return 1;
    60.     }
    61.     ili read(double &x)
    62.     {
    63.         int xx = 0; double f = 1.0, fraction = 1.0; char c = gc(); if (!~c) return EOF;
    64.         while (!isdigit(c)) { if (c == '-') f = -1.0; c = gc(); }
    65.         while (isdigit(c)) { xx = (xx << 3) + (xx << 1) + (c ^ 48), c = gc(); }
    66.         x = xx;
    67.         if (c ^ '.') { x = f * xx; return 1; }
    68.         c = gc();
    69.         while (isdigit(c)) x += (c ^ 48) * (fraction /= 10), c = gc();
    70.         x *= f; return 1;
    71.     }
    72.     ilv write(int x) { if (x < 0) pc('-'), x = -x; if (x > 9) write(x / 10); pc(x % 10 + 48); }
    73.     ilv writeln(int x) { write(x);pc(10); }
    74.     ili read(char *x)
    75.     {
    76.         char c = gc(); if (!~c) return EOF;
    77.         while(!isalpha(c) && !isdigit(c)) c = gc();
    78.         while (isalpha(c) || isdigit(c)) *x++ = c, c = gc();
    79.         *x = 0; return 1;
    80.     }
    81.     ili readln(char *x)
    82.     {
    83.         char c = gc(); if (!~c) return EOF;
    84.         while(c == 10) c = gc();
    85.         while(c >= 32 && c <= 126) *x++ = c, c = gc();
    86.         *x = 0; return 1;
    87.     }
    88.     ilv write(char *x) { while(*x) pc(*x++); }
    89.     ilv write(const char *x) { while(*x) pc(*x++); }
    90.     ilv writeln(char *x) { write(x); pc(10); }
    91.     ilv writeln(const char *x) { write(x); pc(10); }
    92.     ilv write(char c) { pc(c); }
    93.     ilv writeln(char c) { write(c); pc(10); }
    94. } using namespace fastio;
    95. const int maxn = 1e4+5;
    96. int n, m;
    97. struct Point
    98. {
    99.     double x, y;
    100.     Point(double x = 0, double y = 0): x(x), y(y){};
    101.     Point operator - (Point o)
    102.     {
    103.         return Point(x - o.x, y - o.y);
    104.     }
    105.     double operator / (Point o)
    106.     {
    107.         return x * o.y - y * o.x;
    108.     }
    109.     double operator * (Point o)
    110.     {
    111.         return x * o.x + y * o.y;
    112.     }
    113.     Point neg()
    114.     {
    115.         return Point(-x, -y);
    116.     }
    117.     double magnitude()
    118.     {
    119.         return sqrt(SQUARE(x) + SQUARE(y));
    120.     }
    121.     Point scalar(double a)
    122.     {
    123.         return Point(x * a, y *a);
    124.     }
    125.     Point normalize()
    126.     {
    127.         return scalar(1 / magnitude());
    128.     }
    129. } shape1[maxn], shape2[maxn], d;
    130. stack<Point> s;
    131. Point center(Point *shape, int n)
    132. {
    133.     Point ans;
    134.     for (re i = 0; i < n; i++)
    135.     {
    136.         ans.x += shape[i].x;
    137.         ans.y += shape[i].y;
    138.     }
    139.     ans.x /= n, ans.y /= n;
    140.     return ans;
    141. }
    142. Point support1(Point *shape, int n, Point d)
    143. {
    144.     double mx = -inf, proj;
    145.     Point ans;
    146.     for (re i = 0; i < n; i++)
    147.     {
    148.         proj = shape[i] * d;
    149.         if (mx < proj)
    150.         {
    151.             mx = proj;
    152.             ans = shape[i];
    153.         }
    154.     }
    155.     return ans;
    156. }
    157. Point support(Point *shape1, Point *shape2, int n1, int n2, Point d)
    158. {
    159.     Point x = support1(shape1, n1, d), y = support1(shape2, n2, d.neg());
    160.     return x - y;
    161. }
    162. ili dcmp(double x)
    163. {
    164.     if (fabs(x) < eps) return 0;
    165.     return x < 0 ? -1 : 1;
    166. }
    167. Point perp(Point &a, Point &b, Point &c)
    168. {
    169.     return b.scalar(a * c) - a.scalar(b * c);
    170. }
    171. Point closestPointToOrigin(Point &a, Point &b)
    172. {
    173.     double da = a.magnitude();
    174.     double db = b.magnitude();
    175.     double dis = fabs(a / b) / (a - b).magnitude();
    176.     Point ab = b - a, ba = a - b, ao = a.neg(), bo = b.neg();
    177.     if (ab * ao > 0 && ba * bo > 0) return perp(ab, ao, ab).normalize().scalar(dis);
    178.     return da < db ? a.neg() : b.neg();
    179. }
    180. ild gjk(Point *shape1, Point *shape2, int n1, int n2)
    181. {
    182.     d = center(shape2, n2) - center(shape1, n1);
    183.     Point a = support(shape1, shape2, n1, n2, d);
    184.     Point b = support(shape1, shape2, n1, n2, d.neg());
    185.     Point c, p1, p2;
    186.     d = closestPointToOrigin(a, b);
    187.     s.push(a);
    188.     s.push(b);
    189.     while (1)
    190.     {
    191.         c = support(shape1, shape2, n1, n2, d);
    192.         a = s.top(); s.pop();
    193.         b = s.top(); s.pop();
    194.         double da = d * a, db = d * b, dc = d * c;
    195.         if (!dcmp(dc - da) || !dcmp(dc - db)) return d.magnitude();
    196.         p1 = closestPointToOrigin(a, c);
    197.         p2 = closestPointToOrigin(b, c);
    198.         p1.magnitude() < p2.magnitude() ? s.push(a), d = p1 : s.push(b), d = p2;
    199.         s.push(c);
    200.     }
    201. }
    202. signed main()
    203. {
    204. #ifdef LOCAL
    205.     FILE *ALLIN = freopen("d:\\data.in", "r", stdin);
    206. //  freopen("d:\\my.out", "w", stdout);
    207. #endif
    208.     while (read(n), read(m), n + m)
    209.     {
    210.         for (re i = 0; i < n; i++) read(shape1[i].x), read(shape1[i].y);
    211.         for (re i = 0; i < m; i++) read(shape2[i].x), read(shape2[i].y);
    212.         printf("%.5lf\n", gjk(shape1, shape2, n, m));
    213.     }
    214.     flush();
    215. #ifdef LOCAL
    216.     fclose(ALLIN);
    217. #endif
    218.     return 0;
    219. }
    220. https://cloud.tencent.com/developer/article/1679020?areaId=106001
    复制代码

     

     

     

     

    GJK算法计算凸多边形之间的距离
    您需要登录后才可以回帖 登录 | 立即注册

    本版积分规则

    QQ|Archiver|中国膜结构网|中国膜结构协会|进口膜材|国产膜材|ETFE|PVDF|PTFE|设计|施工|安装|车棚|看台|污水池|中国膜结构网_中国空间膜结构协会

    GMT+8, 2024-11-1 08:37 , Processed in 0.151176 second(s), 28 queries .

    Powered by Discuz! X3.5

    © 2001-2024 Discuz! Team.

    快速回复 返回顶部 返回列表