;******************************************************************************************* ;* Double Cusp Group, by Paul Nylander, bugman123.com, 8/6/08 ;******************************************************************************************* ;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 negate (x) (mapcar '- x)) ;Complex Functions (defun complex (x y) (list x y)) (setq I (complex 0.0 1.0)) (defun Re (z) (car z)) (defun Im (z) (cadr z)) (defun conjugate (z) (complex (Re z) (- (Im z)))) (defun CAbs (z) (sqrt (+ (² (Im z)) (² (Re z))))) (defun arg (z) (atan (Im z) (Re z))) (defun CSqr (z) (complex (- (² (Re z)) (² (Im z))) (* 2.0 (Re z) (Im z)))) (defun pow (z n / r theta) (setq r (CAbs z)) (if (= r 0.0) (complex 0.0 0.0) (progn (setq theta (* n (arg z))) (× (complex (cos theta) (sin theta)) (expt r n)) ));i );d (defun CSqrt (z) (pow z 0.5)) (defun × (z1 z2 / n1 n2) (setq n1 (numberp z1) n2 (numberp z2)) (cond ((and n1 n2) (* z1 z2)) (n1 (complex (* z1 (Re z2)) (* z1 (Im z2)))) (n2 (complex (* (Re z1) z2) (* (Im z1) z2))) (T (complex (- (* (Re z1) (Re z2)) (* (Im z1) (Im z2))) (+ (* (Im z1) (Re z2)) (* (Re z1) (Im z2))))) ));d (defun ÷ (z1 z2 / n1 n2) (setq n1 (numberp z1) n2 (numberp z2)) (cond ((and n1 n2) (/ z1 z2)) (n2 (complex (/ (Re z1) z2) (/ (Im z1) z2))) (T (× z1 (pow z2 -1.0))) ));d ;2×2 Complex Matrix Functions (defun matrix (a b c d) (list (list a b) (list c d))) (defun conjugate2 (A) (matrix (conjugate (caar A)) (conjugate (cadar A)) (conjugate (caadr A)) (conjugate (cadadr A)))) (defun det (A) (subtract (× (caar A) (cadadr A)) (× (cadar A) (caadr A)))) (defun mult2 (z A) (matrix (× (caar A) z) (× (cadar A) z) (× (caadr A) z) (× (cadadr A) z))) (defun div2 (A z) (matrix (÷ (caar A) z) (÷ (cadar A) z) (÷ (caadr A) z) (÷ (cadadr A) z))) (defun inverse (A) (div2 (matrix (cadadr A) (negate (cadar A)) (negate (caadr A)) (caar A)) (det A))) (defun dot (A1 A2) (matrix (add (× (caar A1) (caar A2)) (× (cadar A1) (caadr A2))) (add (× (caar A1) (cadar A2)) (× (cadar A1) (cadadr A2))) (add (× (caadr A1) (caar A2)) (× (cadadr A1) (caadr A2))) (add (× (caadr A1) (cadar A2)) (× (cadadr A1) (cadadr A2))) ));d (defun pow2 (A1 n / A2) (setq A2 A1) (repeat (1- n) (setq A2 (dot A2 A1))) A2) (defun CFix (A) (÷ (subtract (subtract (caar A) (cadadr A)) (CSqrt (add (× 4.0 (× (cadar A) (caadr A))) (CSqr (subtract (caar A) (cadadr A)))))) (× 2.0 (caadr A)))) ;Chain Structure (defun chain (z r) (list z r)) (defun zchain (c) (car c)) (defun rchain (c) (cadr c)) (defun ToMatrix (z r) (matrix (÷ z r) (complex (- r (/ (+ (² (Re z)) (² (Im z))) r)) 0.0) (complex (/ 1.0 r) 0.0) (negate (÷ (conjugate z) r)))) ;returns matrix representation of circular chain containing the fix points of M1,M2,M3 (defun MotherCircle (M1 M2 M3 / z0 z1 z2 z3 x1 x2 x3 y1 y2 y3) (setq z1 (CFix M1) x1 (Re z1) y1 (Im z1)) (setq z2 (CFix M2) x2 (Re z2) y2 (Im z2)) (setq z3 (CFix M3) x3 (Re z3) y3 (Im z3)) (setq z0 (÷ (complex (+ (* (² x3) (- y1 y2)) (* (+ (² x1) (* (- y1 y2) (- y1 y3))) (- y2 y3)) (* (² x2) (- y3 y1))) (+ (* -1.0 (² x2) x3) (* (² x1) (- x3 x2)) (* x3 (- y1 y2) (+ y1 y2)) (* x1 (- (+ (² x2) (² y2)) (² x3) (² y3))) (* x2 (- (+ (² x3) (² y3)) (² y1)))) );c (* 2.0 (+ (* x3 (- y1 y2)) (* x1 (- y2 y3)) (* x2 (- y3 y1)))) ));s (ToMatrix z0 (CAbs (subtract z1 z0))) );d ;Calculations (setq ta (complex 1.958591030 -0.011278560) tb (complex 2.0 0.0) tab (÷ (add (× ta tb) (CSqrt (subtract (× (CSqr ta) (CSqr tb)) (× 4.0 (add (CSqr ta) (CSqr tb)))))) 2.0) z0 (÷ (× (subtract tab (complex 2.0 0.0)) tb) (add (subtract (× tb tab) (× 2.0 ta)) (× 2.0 (× I tab)))) b (matrix (÷ (subtract tb (× 2.0 I)) 2.0) (÷ tb 2.0) (÷ tb 2.0) (÷ (add tb (× 2.0 I)) 2.0)) Bi (inverse b) a (dot (matrix tab (÷ (subtract tab (complex 2.0 0.0)) z0) (× (add tab (complex 2.0 0.0)) z0) tab) Bi) Ai (inverse a) c1 (MotherCircle b (dot a (dot b Ai)) (dot a (dot b (dot Ai Bi)))) c2 (MotherCircle (dot b (pow2 a 15)) (dot a (dot b (pow2 a 14))) (dot a (dot b (dot Ai Bi)))) );s ;Create Double Cusp Group (defun reflect (c M) (dot M (dot c (inverse (conjugate2 M))))) (defun zcen (A) (÷ (caar A) (caadr A))) (defun rad (A) (abs (/ 1.0 (Re (caadr A))))) (defun circle (p r color) (entmake (list '(0 . "CIRCLE") (cons 10 p) (cons 40 r) (cons 62 color)))) (defun spiral (c0 M n1 n2 / orbits c) (setq c c0 orbits (list c)) (repeat n1 (setq c (reflect c M)) (if (< (rad c) 100.0) (setq orbits (append (list c) orbits))) );r (setq M (inverse M) c c0) (repeat n2 (setq c (reflect c M)) (if (< (rad c) 100.0) (setq orbits (append orbits (list c)))) );r orbits );d (foreach c (spiral c1 a 83 83) (circle (zcen c) (rad c) 7)) (foreach c (spiral c2 a 91 76) (circle (zcen c) (rad c) 1))