多边形Vonoroi 等分

chlh_jd · 发表于 2012-11-9 00:03

本帖最后由 chlh_jd 于 2012-11-9 00:03 编辑

大家好！
有个几何问题，一直困恼了我很长时间，直到现在都没有解决；

求多边形内的N个点，这N个点构成的Vonoroi形面积为多边形面积的N等分。
如附件所传，求多边形的Vonoroi形7等分的7个点。

这个问题可能有很多答案，只要能求解出一种，其他的方法应该是可以递推的。

chlh_jd · 发表于 2012-11-12 02:36

本帖最后由 chlh_jd 于 2012-11-12 03:14 编辑

写了个近似之一的，希望能有更多几何算法高手参与这个讨论

;;;--------------------------------------------------------------------------------;;;
;;; main function
(defun foo (l n / a cen a1 a2 a3 l1 l2 l3 is p0 ang an1 an d p)
;; middle one , rest distributed in the surrounding method .
;; by GSLS(SS) 2011-11-12
(setq a (f1 l)
cen (car a)
a (cadr a)
a1 (/ a (1- n)))
(setq l1 (list (car l))
l2 (append (cdr l) l1))
;; though Centroid , Aliquot polygon by (1- n) .
;; first point is the startpoint of polygon .
(repeat (- n 2)
(setq l3 nil
is nil
l3 (list (car l1) cen)
a2 (f2 (cons (car l2) l3))
)
(cond ((equal a1 a2 1e-4)
(setq is nil
l1 (cons (car l2) l1)
l2 (cdr l2)))
((> a2 a1)
(setq is T)
)
((< a2 a1)
(while (< a2 a1)
(setq l3 (cons (car l2) l3)
l2 (cdr l2)
a2 (f2 (cons (car l2) l3))))
(if (equal a1 a2 1e-4)
(setq is nil
l1 (cons (car l2) l1)
l2 (cdr l2))
(setq is T))
))
(if is
(progn
(setq p0 (car l2)
ang (angle p0 (car l3))
a3 (- a2 a1)
an1 (angle p0 cen)
an (abs (- ang an1))
)
(if (> an pi)
(setq an (- (+ pi pi) an)))
(setq d (/ (* a3 2.) (sin an) (distance p0 cen))
p (polar p0 ang d)
l1 (cons p l1))))
)
;; cal middle-polygon
;;
(setq l2 (reverse (cons (last l1) l1))
l1 (reverse l1)
a3 (- a1 (/ a n))
l2 (mapcar (function
(lambda (p0 p1 / an)
(setq an (abs (- (angle cen p0) (angle cen p1))))
(if (> an pi)
(setq an (- (+ pi pi) an)))
(/ a3 0.5 (sin an))
)
)
l2
(cdr l2)))
;; solve XY=Z
(setq l2 (solve-xy=z l2 1e-3))
;; list one
(setq l2 (mapcar (function (lambda (p a)
(polar cen (angle cen p) a)))
l1
l2))
(list l1 l2)
)
;;----------
(defun f1 (lst / a c x y a1 c1)
;; get polygon centroid and area
(setq a 0.
c (list 0. 0.)
lst (cons (last lst) lst)
)
(while (cadr lst)
(setq x (car lst)
y (cadr lst)
lst (cdr lst)
a1 (/ (- (* (car x) (cadr y)) (* (car y) (cadr x))) 2.)
a (+ a a1)
c1 (list (/ (* (+ (car x) (car y)) a1) 3.)
(/ (* (+ (cadr x) (cadr y)) a1) 3.))
c (mapcar (function +) c c1)
)
)
(if (= a 0.)
nil
(list (list (/ (car c) a) (/ (cadr c) a)) (abs a))
)
)
;;;-----------
(defun f2 (l)
;;get polygon area
(abs
(/
(apply
(function +)
(mapcar
(function (lambda (x y)
(- (* (car x) (cadr y)) (* (car y) (cadr x)))
))
(cons (last l) l)
l
)
)
2.))
)
;;;
;; (solve-xy=z (list 5 6 7 8) 1e-7)
;; (solve-xy=z '(8.89948e+007 7.06109e+007 7.46468e+007 9.56032e+007 6.8584e+007 7.87004e+007) 1e-7)
(defun solve-xy=z (l eps / f a0 a res)
;; solve Cycle Equations XY=Z , like
;; a*b = A ; b*c = B ; c*d = C ; d*e = D ; e*f = E ; f*a = F .
;; Only for polygon Vonoroi-devide
;; by GSLS(SS) 2012-11-12
(defun f (a0 l)
(foreach b l
(setq a0 (/ b a0))
)
)
(setq a0 (sqrt (car l))
a (f a0 l))
(while (not (equal a0 a eps))
(setq a0 (sqrt (* a0 a))
a (f a0 l))
(if (or (> a0 (* (sqrt (car l)) 3.16227766))
(< a0 (/ (sqrt (car l)) 3.16227766)))
(setq a0 (sqrt (car l))
a (f a0 l)
eps (* 10. eps))))
(setq res (list (sqrt (* a a0))))
(foreach b l
(setq a (/ b (car res))
res (cons a res)))
(reverse (cdr res))
)
;;-------------
(defun ss-assoc (a lst / b res)
(while (setq b (assoc a lst))
(setq lst (cdr (member b lst))
res (cons (cdr b) res)
))(reverse res))
;;
;;--------------------------------
;; for test
(defun c:test (/ pl ent l n res)
(if (and
(setq pl (car (entsel "\nSelect a LWPolyLine :")))
(setq ent (entget pl)
l (ss-assoc 10 ent))
(setq n (getint "\n Vonoroi devide Number : ")))
(if (setq res (foo l n))
(progn
(setq l1 (car res)
l2 (cadr res))
(mapcar (function (lambda (p1 p2)
(entmake (list (cons 0 "LINE")
(cons 10 p1)
(cons 11 p2)
(cons 62 1)
(cons 8 "Temp")))))
(append l1 l2)
(append l2 (cdr l2) (list (car l2))))
)
)
)
(princ)
)

;;;--------------------------------------------------------------------------------;;; 
;;; main function 
(defun foo  (l n / a cen a1 a2 a3 l1 l2 l3 is p0 ang an1 an d p) 
  ;; middle one ,  rest distributed in the surrounding method . 
  ;; by GSLS(SS) 2011-11-12 
  (setq a   (f1 l) 
 cen (car a) 
 a   (cadr a) 
 a1  (/ a (1- n))) 
  (setq l1 (list (car l)) 
 l2 (append (cdr l) l1)) 
  ;; though Centroid , Aliquot polygon by (1- n) . 
  ;; first point is the startpoint of polygon . 
  (repeat (- n 2) 
    (setq l3 nil 
   is nil 
   l3 (list (car l1) cen) 
   a2 (f2 (cons (car l2) l3)) 
   ) 
    (cond ((equal a1 a2 1e-4) 
    (setq is nil 
   l1 (cons (car l2) l1) 
   l2 (cdr l2))) 
   ((> a2 a1) 
    (setq is T) 
    ) 
   ((< a2 a1) 
    (while (< a2 a1) 
      (setq l3 (cons (car l2) l3) 
     l2 (cdr l2) 
     a2 (f2 (cons (car l2) l3)))) 
    (if (equal a1 a2 1e-4) 
      (setq is nil 
     l1 (cons (car l2) l1) 
     l2 (cdr l2)) 
      (setq is T)) 
    )) 
    (if is 
      (progn 
 (setq p0  (car l2) 
       ang (angle p0 (car l3)) 
       a3  (- a2 a1) 
       an1 (angle p0 cen) 
       an  (abs (- ang an1)) 
       ) 
 (if (> an pi) 
   (setq an (- (+ pi pi) an))) 
 (setq d  (/ (* a3 2.) (sin an) (distance p0 cen)) 
       p  (polar p0 ang d) 
       l1 (cons p l1)))) 
    ) 
  ;; cal middle-polygon 
  ;; 
  (setq l2 (reverse (cons (last l1) l1)) 
 l1 (reverse l1) 
 a3 (- a1 (/ a n)) 
 l2 (mapcar (function 
       (lambda (p0 p1 / an) 
         (setq an (abs (- (angle cen p0) (angle cen p1)))) 
         (if (> an pi) 
    (setq an (- (+ pi pi) an))) 
         (/ a3 0.5 (sin an)) 
         ) 
       ) 
     l2 
     (cdr l2))) 
  ;; solve XY=Z 
  (setq l2 (solve-xy=z l2 1e-3)) 
  ;; list one 
  (setq l2 (mapcar (function (lambda (p a) 
          (polar cen (angle cen p) a))) 
     l1 
     l2)) 
  (list l1 l2) 
  ) 
;;---------- 
(defun f1  (lst / a c x y a1 c1) 
  ;; get polygon centroid and area 
  (setq a   0. 
 c   (list 0. 0.) 
 lst (cons (last lst) lst) 
 ) 
  (while (cadr lst) 
    (setq x   (car lst) 
   y   (cadr lst) 
   lst (cdr lst) 
   a1  (/ (- (* (car x) (cadr y)) (* (car y) (cadr x))) 2.) 
   a   (+ a a1) 
   c1  (list (/ (* (+ (car x) (car y)) a1) 3.) 
      (/ (* (+ (cadr x) (cadr y)) a1) 3.)) 
   c   (mapcar (function +) c c1) 
   ) 
    ) 
  (if (= a 0.) 
    nil 
    (list (list (/ (car c) a) (/ (cadr c) a)) (abs a)) 
    ) 
  ) 
;;;----------- 
(defun f2  (l) 
  ;;get polygon area 
  (abs 
    (/ 
      (apply 
 (function +) 
 (mapcar 
   (function (lambda (x y) 
        (- (* (car x) (cadr y)) (* (car y) (cadr x))) 
        )) 
   (cons (last l) l) 
   l 
   ) 
 ) 
      2.)) 
  ) 
;;; 
;; (solve-xy=z (list 5 6 7 8) 1e-7) 
;; (solve-xy=z '(8.89948e+007  7.06109e+007 7.46468e+007 9.56032e+007 6.8584e+007 7.87004e+007) 1e-7) 
(defun solve-xy=z  (l eps / f  a0 a res) 
  ;; solve Cycle Equations XY=Z , like 
  ;;  a*b = A ; b*c = B ; c*d = C ; d*e = D ; e*f = E ; f*a = F . 
  ;; Only for polygon Vonoroi-devide 
  ;; by GSLS(SS) 2012-11-12 
  (defun f  (a0 l) 
    (foreach b l 
      (setq a0 (/ b a0)) 
      ) 
    ) 
  (setq a0 (sqrt (car l)) 
 a  (f a0 l)) 
  (while (not (equal a0 a eps)) 
    (setq a0 (sqrt (* a0 a)) 
   a  (f a0 l)) 
    (if (or (> a0 (* (sqrt (car l)) 3.16227766)) 
     (< a0 (/ (sqrt (car l)) 3.16227766))) 
      (setq a0 (sqrt (car l)) 
     a (f a0 l) 
     eps (* 10. eps)))) 
  (setq res (list (sqrt (* a a0)))) 
  (foreach b  l 
    (setq a   (/ b (car res)) 
   res (cons a res))) 
(reverse (cdr res))    
  ) 
;;------------- 
(defun ss-assoc (a lst / b res) 
  (while (setq b (assoc a lst)) 
    (setq lst  (cdr (member b lst)) 
   res (cons (cdr b) res) 
    ))(reverse res)) 
;; 
;;-------------------------------- 
;; for test 
(defun c:test  (/ pl ent l n res) 
  (if (and 
 (setq pl (car (entsel "\nSelect a LWPolyLine :"))) 
 (setq ent (entget pl) 
       l   (ss-assoc 10 ent)) 
 (setq n (getint "\n Vonoroi devide Number : "))) 
    (if (setq res (foo l n)) 
      (progn 
 (setq l1 (car res) 
       l2 (cadr res)) 
 (mapcar (function (lambda (p1 p2) 
       (entmake (list (cons 0 "LINE") 
        (cons 10 p1) 
        (cons 11 p2) 
        (cons 62 1) 
        (cons 8 "Temp"))))) 
  (append l1 l2) 
  (append l2 (cdr l2) (list (car l2)))) 
 ) 
      ) 
    ) 
  (princ) 
  )

qjchen · 发表于 2012-11-15 23:04

前些日子看了，不过不大懂

也没有查到太类似的网页

唯一的一个是 rhino的grasshopper的一个插件，似乎是提供了源码，但似乎只是一个描述语言，看不懂

http://crtl-i.com/blog/2010/06/same-area-voronoi-using-galapagos/

另外，查到的类似文章，人家经常做的是，已有一些点，然后，需要做一些图形划分，每个图形中包含一个点，且面积相等

比但高山流水朋友的是要voronoi，点是不确定的，说不定更难，一时也想不到什么好方法

chlh_jd · 发表于 2012-11-18 00:13

本帖最后由 chlh_jd 于 2012-11-18 00:27 编辑

非常感谢QJChen老师！
每次总能得到您的帮助

曲面有限元分析程序经常会用到自适应网格划分，如果是自适应Delaunay三角网划分，实际上就得到了Vonoroi 图，但是麻烦的是这个自适应划分和面积相等可能对于任何一个可行方案都得计算很长时间。

国外一些先进的复杂结构模型建模辅助工具（如“犀牛”之类的参数化建模）可能有些帮助，但我一直没能得到这些软件试用

myqzq · 发表于 2013-4-27 19:52

厉害，很强大，收藏了

freehand8008 · 发表于 2013-6-28 10:02

chlh_jd 发表于 2012-11-12 02:36
写了个近似之一的，希望能有更多几何算法高手参与这个讨论

defun foo ?
输入 foo 为什么不行啊

434939575 · 发表于 2014-12-30 11:40

高级活，长见识了！谢谢分享！

看天的小树 · 发表于 2015-6-3 19:59

freehand8008 发表于 2013-6-28 10:02
defun foo ?
输入 foo 为什么不行啊

你没看清楚前面没有带 C：啊

只能外部调用它

看天的小树 · 发表于 2015-6-3 20:02

freehand8008 发表于 2013-6-28 10:02
defun foo ?
输入 foo 为什么不行啊

输入test，那个是主函数

freehand8008 · 发表于 2015-6-4 09:34

看天的小树发表于 2015-6-3 20:02
输入test，那个是主函数

成功

账号		自动登录	找回密码
密码			注册

多边形Vonoroi 等分

本帖子中包含更多资源

本帖被以下淘专辑推荐:

本帖子中包含更多资源

点评

评分

本帖子中包含更多资源