明经CAD社区

 找回密码
 注册

QQ登录

只需一步,快速开始

搜索
查看: 3276|回复: 7

[自我挑战] 用鲍威尔法求解最优化问题

  [复制链接]
发表于 2013-10-6 16:12:24 | 显示全部楼层 |阅读模式
本帖最后由 chlh_jd 于 2013-10-8 17:44 编辑

大家好!
很久以前就在想:用LISP写一个求解最优化的通用函数,因为LISP语言本身就比较适合;基于对数学的弱弱了解,写了个粗糙的函数和大家讨论,希望您能多参与这个话题,感谢在先!

  1. ;;(Powell val con arg arg_limit eps)
  2. (defun Powell  (val                  con                    arg
  3.                 arg_limit          eps                    / ss_powell_arg0
  4.                 ss_powell_1dmin          ss_powell_main_function
  5.                 ss_powell_a0          ss_powell_dfm_lst ss_powell_dfm
  6.                 ss_powell_x0          ss_powell_f0            ss_powell_x1
  7.                 ss_powell_f1          ss_powell_x2            ss_powell_f2)
  8.   ;; --------------------------------------------------------;;
  9.   ;; Use Powell method to solve solve optimization problems
  10.   ;;
  11.   ;; Args:
  12.   ;; val -- Expressions calculation function
  13.   ;; con -- Restriction function
  14.   ;; arg -- A string list of the Parameters's name , such as '( "x1" "x2" ...)
  15.   ;; arg_limit -- parameters's limits , like '( (-1e14 1e14) (0 1e14) ...)
  16.   ;; eps -- accuracy
  17.   ;;
  18.   ;; by GSLS(SS) 2013.10.6
  19.   ;;---------------------------------------------------------;;
  20.   (defun ss_powell_1dmin  (x -xm x0 +xm / f0 fa fb -fm +fm)
  21.     ;; One-dimensional search
  22.     ;; This function has not the Powell method .
  23.     ;; How to change it connect with Gradient direction ?
  24.     (set (read x) x0)
  25.     (setq f0 (if (con)
  26.                (val)
  27.                1e99))
  28.     (set (read x) -xm)
  29.     (setq fa (if (con)
  30.                (val)
  31.                1e99))
  32.     (set (read x) (* 0.5 (+ x0 -xm)))
  33.     (setq -fm (if (con)
  34.                 (val)
  35.                 1e99))
  36.     (set (read x) +xm)
  37.     (setq fb (if (con)
  38.                (val)
  39.                1e99))
  40.     (set (read x) (* 0.5 (+ x0 +xm)))
  41.     (setq +fm (if (con)
  42.                 (val)
  43.                 1e99))
  44.     (cond
  45.       ((or (equal x0 -xm eps) (equal x0 +xm eps))
  46.        x0)
  47.       ((< fa -fm f0 +fm fb)
  48.        -xm)
  49.       ((> fa -fm f0 +fm fb)
  50.        +xm)
  51.       ((and (> fa -fm f0) (< f0 +fm fb))
  52.        (ss_powell_1dmin x (*  0.5 (+ x0 -xm)) x0 (* 0.5 (+ x0 +xm))))
  53.       ((and (> fa -fm) (< -fm f0 +fm fb))
  54.        (ss_powell_1dmin x -xm (* 0.5 (+ x0 -xm)) x0))
  55.       ((and (> fb +fm) (> fa -fm f0 +fm))
  56.        (ss_powell_1dmin x x0 (* 0.5 (+ x0 +xm)) +xm))
  57.       ((and (equal x0 -xm (sqrt eps)) (equal x0 +xm (sqrt eps)))
  58.        (/ (+ x0 -xm +xm) 3.))
  59.       ))
  60.   (setq ss_powell_arg0 arg)
  61.   (if arg_limit
  62.     (mapcar
  63.       (function
  64.         (lambda        (x y)
  65.           (cond
  66.             ((or (not y) (apply (function or) y))
  67.              (set (read (strcat x "_limit")) (list -1e14 1e14)))
  68.             ((not (car y))
  69.              (set (read (strcat x "_limit")) (list -1e14 (cadr y))))
  70.             ((not (cadr y))
  71.              (set (read (strcat x "_limit")) (list (car y) 1e14)))
  72.             (t (set (read (strcat x "_limit")) y))
  73.             )))
  74.       arg
  75.       arg_limit)
  76.     (foreach x        arg
  77.       (set (read (strcat x "_limit")) (list -1e14 1e14)))
  78.     ) ;_ x1_limit
  79.   (setq        ss_powell_x0
  80.         (mapcar
  81.            (function
  82.              (lambda (x)
  83.                (*
  84.                 0.5
  85.                 (+
  86.                    (car (eval (read (strcat x "_limit"))))
  87.                    (cadr (eval (read (strcat x "_limit"))))))))
  88.            arg))  
  89.   (defun ss_powell_main_function  (/ ss_powell_i)
  90.     (ss-set arg ss_powell_x0)
  91.     (setq ss_powell_f0
  92.            (if (con)
  93.              (val)
  94.              1e99))
  95.     (setq ss_powell_a0
  96.            (mapcar (function (lambda (x)
  97.                                0.0))
  98.                    arg)
  99.           ss_powell_dfm_lst
  100.            nil
  101.           ss_powell_i 0)
  102.     (mapcar
  103.       (function
  104.         (lambda        (x  /
  105.                 ss_powell_a1_f0      ss_powell_a1_x0
  106.                 ss_powell_a1_l              ss_powell_a1_-xm
  107.                 ss_powell_a1_+xm     ss_powell_a1_ai
  108.                 ss_powell_a1_f1      )
  109.           (setq        ss_powell_a1_f0
  110.                 (if(con)(val)1e99))
  111.           (setq ss_powell_a1_x0 (eval (read x)))
  112.           (setq ss_powell_a1_l (eval (read (strcat x "_limit"))))
  113.           (setq        ss_powell_a1_-xm (car ss_powell_a1_l)
  114.                 ss_powell_a1_+xm (cadr ss_powell_a1_l))
  115.           (setq        ss_powell_a1_ai
  116.                 (- (ss_powell_1dmin
  117.                       x
  118.                       ss_powell_a1_-xm
  119.                       ss_powell_a1_x0
  120.                       ss_powell_a1_+xm)
  121.                     ss_powell_a1_x0))
  122.           (setq ss_powell_a0 (ch-lst ss_powell_a1_ai ss_powell_i ss_powell_a0))
  123.           (ss-set arg
  124.                   (mapcar (function +) ss_powell_x0 ss_powell_a0))
  125.           (setq        ss_powell_a1_f1
  126.                 (if (con)
  127.                    (val)
  128.                    1e99))         
  129.           (setq        ss_powell_dfm_lst
  130.                 (cons (- ss_powell_a1_f0 ss_powell_a1_f1)        ss_powell_dfm_lst))
  131.           (setq ss_powell_i (1+ ss_powell_i))         
  132.           ))
  133.       arg)
  134.     (setq ss_powell_dfm_lst (reverse ss_powell_dfm_lst))
  135.     (setq ss_powell_x1
  136.            (mapcar (function +) ss_powell_x0 ss_powell_a0))
  137.     (setq ss_powell_x2
  138.            (mapcar (function (lambda (x y)
  139.                                (- (* 2. x) y)))
  140.                    ss_powell_x1
  141.                    ss_powell_x0))
  142.     (ss-set arg ss_powell_x1)
  143.     (setq ss_powell_f1
  144.            (if (con)
  145.              (val)
  146.              1e99))
  147.     (ss-set arg ss_powell_x2)
  148.     (setq ss_powell_f2
  149.            (if (con)
  150.              (val)
  151.              1e99)))  
  152.   (ss_powell_main_function)
  153.   (while (not (equal ss_powell_x1 ss_powell_x0 eps))
  154.     (setq ss_powell_dfm (apply (function max) ss_powell_dfm_lst))
  155.     (if
  156.       (or (>= ss_powell_f2 ss_powell_f0)
  157.           (>= (* (+ ss_powell_f0 ss_powell_f2 (* -2. ss_powell_f1))
  158.                 (- ss_powell_f0 ss_powell_f1 ss_powell_dfm))
  159.               (* 0.5
  160.                 ss_powell_dfm
  161.                 (expt (- ss_powell_f0 ss_powell_f2) 2))))
  162.       (progn
  163.        (if (< ss_powell_f1 ss_powell_f2)
  164.         (setq ss_powell_x0 ss_powell_x1)
  165.         (setq ss_powell_x0 ss_powell_x2))
  166.       ;_ (setq ss_powell_dfm_i (vl-sort-i ss_powell_dfm_lst (function >)))   
  167.     )
  168.       (progn;_this is not
  169.         (setq ss_powell_i
  170.                 (car
  171.                   (vl-sort-i ss_powell_dfm_lst (function >))))        
  172.         (setq ss_powell_a0
  173.                 ((lambda(/ i)
  174.                   (setq i -1)
  175.                   (mapcar
  176.                     (function (lambda (y)
  177.                                 (setq i (1+ i))
  178.                                 (if (= i ss_powell_i)
  179.                                   y
  180.                                   (* -2.  y))))                    
  181.                     ss_powell_a0))))
  182.         (setq ss_powell_x0 (mapcar (function +) ss_powell_x0 ss_powell_a0))
  183.         ;_(setq ss_powell_dfm_i (shift(vl-sort-i ss_powell_dfm_lst (function >))))
  184.         )
  185.       ;_(setq ss_powell_x0 ss_powell_x1)
  186.       )
  187.     ;|
  188.     (setq arg (mapcar (function (lambda (i)
  189.                                   (nth i arg)))
  190.                       ss_powell_dfm_i))
  191.     (setq ss_powell_x0 (mapcar (function (lambda (i)
  192.                                   (nth i ss_powell_x0)))
  193.                       ss_powell_dfm_i))|;
  194.     (ss_powell_main_function))
  195.   (ss-set arg ss_powell_x1)
  196.   (setq        ss_powell_f1
  197.         (if (con)
  198.            (val)
  199.            1e99))
  200.   (list (mapcar (function (lambda (x) (eval(read x)))) ss_powell_arg0) ss_powell_f1))

  201. ;;; Used function
  202. (defun ss-set  (ss_set_str_lst ss_set_lst / ss_set_i)
  203.   (setq ss_set_i 0)
  204.   (repeat (length ss_set_lst)
  205.     (set (read (nth ss_set_i ss_set_str_lst))
  206.   (nth ss_set_i ss_set_lst))
  207.     (setq ss_set_i (1+ ss_set_i))))
  208. ;;-----------------------------
  209. (defun ch-lst  (new i lst / j len fst mid)
  210.   (if (/= (type i) (quote list))
  211.     (cond
  212.       ((not (listp lst))
  213.        lst)
  214.       ((minusp i) lst)
  215.       ((> i (setq len (length lst))) lst)
  216.       ((> i (/ len 2))
  217.        (reverse (ch-lst new (1- (- len i)) (reverse lst))))
  218.       (t
  219.        (append
  220.   (progn
  221.      (setq fst nil)
  222.      (repeat (rem i 4)
  223.        (setq fst (cons (car lst) fst)
  224.        lst (cdr lst)))
  225.      (repeat (/ i 4)
  226.        (setq fst (vl-list* (cadddr lst)
  227.         (caddr lst)
  228.         (cadr lst)
  229.         (car lst)
  230.         fst)
  231.        lst (cddddr lst)))
  232.      (reverse fst))
  233.   (list new)
  234.   (cdr lst))))
  235.     (progn
  236.       (setq j (cadr i)
  237.       i (car i))
  238.       (if j
  239.   (progn
  240.     (setq  mid (nth i lst)
  241.     mid (ch-lst new j mid))
  242.     (ch-lst mid i lst))
  243.   (ch-lst new i lst)))))
应用例子1:求10(x1+x2-5)^2+(x1-x2)^2的极小点和极小值,理论解: ((2.5 2.5) 0)
  1. (defun val  ()
  2.   ;; 10(x1+x2-5)^2+(x1-x2)^2
  3.   (+ (* 10. (expt (+ x1 x2 -5.) 2)) (expt (- x1 x2) 2)))
  4. (defun con () T)
  5. (setq arg  (list "x1" "x2")
  6.       arg_limit  nil
  7.       eps  1e-6)
  8. (Powell val con arg arg_limit eps);;->((2.5 2.5) 1.32989e-011)
应用例子2:求60-10x1-4x2+x1^2+x2^2-x1x2的极小点和极小值,理论解:((8 6) 8)
  1. (defun val  ()
  2.   ;;60-10x1-4x2+x1^2+x2^2-x1x2
  3.   (+ 60. (* -10. x1) (* -4. x2) (* x1 x1) (* x2 x2) (* -1. x1 x2)))
  4. (defun con () T)
  5. (setq arg  (list "x1" "x2")
  6.       arg_limit  nil
  7.       eps  1e-6)
  8. (Powell val con arg arg_limit eps);;->((8.0 6.0) 8.0)

评分

参与人数 1明经币 +1 金钱 +30 收起 理由
highflybird + 1 + 30

查看全部评分

发表于 2013-10-8 08:34:02 | 显示全部楼层
这么深的东西,实际工作中从来没有用到过,今后可能也用不到。帮你顶一下
发表于 2013-10-8 12:56:24 | 显示全部楼层
一个字
伪原码

点评

何出此言?  发表于 2022-5-30 12:05
 楼主| 发表于 2013-10-8 17:36:33 | 显示全部楼层
哦,编辑的时候,前面的被我删除了
发表于 2017-11-11 10:41:09 | 显示全部楼层
(defun val  ()
  ;;x1^2+2x2^2-2x1x2-4x1
  (+ (* x1 x1) (* 2. x2 x2) (* -4. x1)  (* -2. x1 x2)))
(defun con () T)
(setq arg  (list "x1" "x2")
      arg_limit  nil
      eps  1e-6)
(Powell val con arg arg_limit eps)

为什么以上运行不出来呢?
发表于 2022-5-27 14:42:25 | 显示全部楼层
Powell用lisp实现很特别。赞一个

发表于 2022-5-27 15:09:01 | 显示全部楼层
楼主函数可以运行。  改变函数后迭代 有时导致acad死机

主程序迭代是否增加约束条件,,否则收敛慢的函数计算时间过长,容易死机
您需要登录后才可以回帖 登录 | 注册

本版积分规则

小黑屋|手机版|CAD论坛|CAD教程|CAD下载|联系我们|关于明经|明经通道 ( 粤ICP备05003914号 )  
©2000-2023 明经通道 版权所有 本站代码,在未取得本站及作者授权的情况下,不得用于商业用途

GMT+8, 2024-11-25 03:35 , Processed in 0.165765 second(s), 26 queries , Gzip On.

Powered by Discuz! X3.4

Copyright © 2001-2021, Tencent Cloud.

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