;******************************************************************************************* ;* Boy's Surface, by Paul Nylander, bugman123.com, 3/11/09 ;******************************************************************************************* ;Basic Functions (defun ² (x) (* x x)) (defun add (x1 x2) (if (listp x1) (mapcar 'add x1 x2) (+ x1 x2))) (defun subtract (x1 x2) (if (listp x1) (mapcar 'subtract x1 x2) (- x1 x2))) (defun ÷ (x1 x2) (if (listp x1) (mapcar '(lambda (x) (÷ x x2)) x1) (/ x1 x2))) (defun norm (x) (apply '+ (mapcar '² x))) ;Complex Functions (defun complex (x y) (list x y)) (setq CR (complex 1.0 0.0) CI (complex 0.0 1.0)) (defun Re (z) (car z)); (if (numberp z) z (car z)) (defun Im (z) (cadr z)); (if (numberp z) 0.0 (cadr z)) (defun CAbs (z) (sqrt (+ (² (Im z)) (² (Re z))))) (defun arg (z) (atan (Im z) (Re z))) (defun CPow (z n / r theta) (setq r (CAbs z)) (if (= r 0.0) (complex 0.0 0.0) (progn (setq theta (* n (arg z))) (CMult (expt r n) (complex (cos theta) (sin theta))) ));i );d (defun CMult (z1 z2) (cond ((numberp z1) (complex (* z1 (Re z2)) (* z1 (Im z2)))) ((numberp z2) (complex (* z2 (Re z1)) (* z2 (Im z1)))) (T (complex (- (* (Re z1) (Re z2)) (* (Im z1) (Im z2))) (+ (* (Im z1) (Re z2)) (* (Re z1) (Im z2))))) ));d (defun CDiv (z1 z2) (CMult z1 (CPow z2 -1.0))) (defun CExp (z) (CMult (exp (Re z)) (complex (cos (Im z)) (sin (Im z))))) ;Converts Polygon into AutoCAD's Triangle Elements, only works for convex polygons ;(TriangulateElem '(1 2 3 4 5 6 7 8 9 10)) -> '((1 2 -3) (-1 3 -4) (-1 4 -5) (-1 5 -6) (-1 6 -7) (-1 7 -8) (-1 8 -9) (-1 9 10)) (defun TriangulateElem (elem / n i i0 i1 i2 elems) (setq n (length elem)) (if (= n 3) (list elem) (progn (setq i0 (car elem) i1 (caddr elem) elems (list (list i0 (cadr elem) (- i1))) i 3) (while (< i (1- n)) (setq i2 (nth i elem) elems (append elems (list (list (- i0) i1 (- i2)))) i1 i2 i (1+ i))) (append elems (list (list (- i0) i1 (last elem)))) ));i );d (defun mesh (nodes elems1 / elems2 nnode nelem elem) (setq elems2 (apply 'append (mapcar 'TriangulateElem elems1))); introduces AutoCAD's negative index notation (setq nnode (length nodes) nelem (length elems2)) (entmake (list '(0 . "POLYLINE") '(70 . 64) (cons 71 nnode) (cons 72 nelem))) (foreach p nodes (entmake (list '(0 . "VERTEX") (cons 10 p) '(70 . 192)))) (foreach elem elems2 (entmake (list '(0 . "VERTEX") '(10 0.0 0.0 0.0) '(70 . 128) (cons 71 (car elem)) (cons 72 (cadr elem)) (cons 73 (caddr elem)) (cons 74 (- (abs (car elem)))) )));f (entmake '((0 . "SEQEND"))) );d (defun ParametricPlot3D (f u1 u2 du v1 v2 dv / u v i j i1 i2 i3 i4 imax jmax nodes elems) (setq imax (fix (+ (/ (- u2 u1) du) 0.5)) jmax (fix (+ (/ (- v2 v1) dv) 0.5))) (setq j 0); calculate nodes (while (<= j jmax) (setq v (+ v1 (* j dv)) i 0) (while (<= i imax) (setq u (+ u1 (* i du)) nodes (append nodes (list (f u v))) i (1+ i))) (setq j (1+ j)) );w (setq i 1); assemble quad elements (while (<= i imax) (setq j 0) (while (< j jmax) (setq i1 (+ i (* j (1+ imax))) i2 (1+ i1) i3 (+ i1 (1+ imax) 1) i4 (+ i1 (1+ imax)) elems (append elems (list (list i1 i2 i3 i4))) j (1+ j) ));w (setq i (1+ i)) );w (mesh nodes elems) );d ;Boy's Surface (defun BoysSurface (r theta / z a m) (setq z (CMult r (CExp (CMult CI theta)))) (setq a (subtract (add (CPow z 6.0) (CMult (sqrt 5.0) (CPow z 3.0))) CR)) (setq m (list (Im (CDiv (subtract (CPow z 5.0) z) a)) (Re (CDiv (add (CPow z 5.0) z) a)) (+ (Im (CDiv (add (CPow z 6.0) CR) (CMult 1.5 a))) 0.5) ));s (mapcar '(lambda (x) (if (= x 0.0) 0.0 x)) (÷ m (norm m))) );d (ParametricPlot3D BoysSurface 1.0e-6 1.0 (/ (- 1.0 1.0e-6) 20) (- pi) pi (/ pi 45)) (command "vpoint" (list -1.0 -1.0 1.0)); (command "zoom" "extents")