Skip to content

Instantly share code, notes, and snippets.

Embed
What would you like to do?
クライン群の極限集合を描画する
(defparameter *inf* (* most-positive-fixnum #c(1 1)))
(defclass matrix ()
((a :accessor mat-a :initform 1 :initarg :a)
(b :accessor mat-b :initform 0 :initarg :b)
(c :accessor mat-c :initform 0 :initarg :c)
(d :accessor mat-d :initform 1 :initarg :d)))
(defgeneric mat* (a b))
(defgeneric inverse (m))
(defgeneric mat-trace (m))
(defmethod mat* ((m1 matrix) (m2 matrix))
(with-accessors ((m1-a mat-a) (m1-b mat-b) (m1-c mat-c) (m1-d mat-d)) m1
(with-accessors ((m2-a mat-a) (m2-b mat-b) (m2-c mat-c) (m2-d mat-d)) m2
(make-instance 'matrix
:a (+ (* m1-a m2-a) (* m1-b m2-c))
:b (+ (* m1-a m2-b) (* m1-b m2-d))
:c (+ (* m1-c m2-a) (* m1-d m2-c))
:d (+ (* m1-c m2-b) (* m1-d m2-d))))))
(defmethod mat* ((coefficient complex) (m1 matrix))
(with-accessors ((m1-a mat-a) (m1-b mat-b) (m1-c mat-c) (m1-d mat-d)) m1
(make-instance 'matrix
:a (* coefficient m1-a)
:b (* coefficient m1-b)
:c (* coefficient m1-c)
:d (* coefficient m1-d))))
(defmethod mat* ((coefficient integer) (m1 matrix))
(with-accessors ((m1-a mat-a) (m1-b mat-b) (m1-c mat-c) (m1-d mat-d)) m1
(make-instance 'matrix
:a (* coefficient m1-a)
:b (* coefficient m1-b)
:c (* coefficient m1-c)
:d (* coefficient m1-d))))
(defmethod inverse ((m matrix))
(let ((determinant (- (* (mat-a m) (mat-d m))
(* (mat-b m) (mat-c m)))))
(if (= determinant 0)
(make-instance 'matrix
:a *inf*
:b *inf*
:c *inf*
:d *inf*)
(make-instance 'matrix
:a (/ (mat-d m) determinant)
:b (/ (* -1 (mat-b m)) determinant)
:c (/ (* -1 (mat-c m)) determinant)
:d (/ (mat-a m) determinant)))))
(defmethod mat-trace ((m matrix))
(+ (mat-a m) (mat-d m)))
(defmethod print-object ((m matrix) s)
(write-char #\{ s)
(princ (mat-a m) s)
(write-char #\, s)
(princ (mat-b m) s)
(write-char #\LineFeed s)
(write-char #\space s)
(princ (mat-c m) s)
(write-char #\, s)
(princ (mat-d m) s)
(write-char #\} s))
(defun fix-plus ((m matrix))
(with-accessors ((a mat-a) (b mat-b) (c mat-c) (d mat-d)) m
(/ (+ (- a d) (sqrt (- (expt (mat-trace m) 2) 4)))
(* c 2))))
(defun mobius-on-point ((m matrix) p)
(with-accessors ((a mat-a) (b mat-b) (c mat-c) (d mat-d)) m
(cond ((= *inf* c) (cond ((= c 0) *inf*)
(t (/ a c))))
(t (let ((numerix (+ (* a p) b))
(denominator (+ (* c p) d)))
(cond ((= denominator 0) *inf*)
(t (/ numerix denominator))))))))
(defun init (max-level epsilon)
(defparameter *fix-point* (make-array '(5 4)))
(defparameter *gens* (make-array '(5)))
(defparameter *max-level* max-level)
(defparameter *word* (make-array (list (+ *max-level* 1))))
(defparameter *tags* (make-array (list (+ *max-level* 1))))
(defparameter *epsilon* epsilon)
(defparameter *line-points-list* ())
(defparameter *level* 1))
(defun print-tags ()
(loop for i from 1 to *level* do
(princ (aref *tags* i)))
(print ""))
(defun set-gens ((t_a complex) (t_b complex))
(let* ((t_ab (/ (- (* t_a t_b)
(sqrt (- (* (expt t_a 2) (expt t_b 2))
(* 4 (+ (expt t_a 2) (expt t_b 2))))))
2))
(z0 (/ (* (- t_ab 2) t_b)
(+ (- (* t_b t_ab) (* t_a 2))
(* t_ab #c(0 2)))))
(a (make-instance 'matrix
:a (/ t_a 2)
:b (/ (+ (- (* t_a t_ab) (* t_b 2)) #c(0 4))
(* z0 (+ (* t_ab 2) 4)))
:c (/ (* z0
(- (* t_a t_ab) (* t_b 2) #c(0 4)))
(- (* t_ab 2) 4))
:d (/ t_a 2)))
(b (make-instance 'matrix
:a (/ (- t_b #c(0 2)) 2)
:b (/ t_b 2)
:c (/ t_b 2)
:d (/ (+ t_b #c(0 2)) 2)))
(_A (inverse a))
(_B (inverse b)))
(setf (aref *fix-point* 1 1) (fix-plus (reduce #'mat* (list b _A _B a))))
(setf (aref *fix-point* 1 2) (fix-plus a))
(setf (aref *fix-point* 1 3) (fix-plus (reduce #'mat* (list _B _A b a))))
(setf (aref *fix-point* 2 1) (fix-plus (reduce #'mat* (list _A _B a b))))
(setf (aref *fix-point* 2 2) (fix-plus b))
(setf (aref *fix-point* 2 3) (fix-plus (reduce #'mat* (list a _B _A b))))
(setf (aref *fix-point* 3 1) (fix-plus (reduce #'mat* (list _B a b _A))))
(setf (aref *fix-point* 3 2) (fix-plus _A))
(setf (aref *fix-point* 3 3) (fix-plus (reduce #'mat* (list b a _B _A))))
(setf (aref *fix-point* 4 1) (fix-plus (reduce #'mat* (list a b _A _B))))
(setf (aref *fix-point* 4 2) (fix-plus _B))
(setf (aref *fix-point* 4 3) (fix-plus (reduce #'mat* (list _A b a _B))))
(setf (aref *gens* 1) a)
(setf (aref *gens* 2) b)
(setf (aref *gens* 3) _A)
(setf (aref *gens* 4) _B)
(setf (aref *word* 0) (make-instance 'matrix))
(setf (aref *word* 1) (aref *gens* 1))
(setf (aref *tags* 1) 1)))
(defun go-forward ()
(if (eq (branch-termination-p) nil)
(progn
(incf *level*)
(setf (aref *tags* *level*)
(abs (mod (+ (aref *tags* (- *level* 1)) 1) 4)))
(if (= (aref *tags* *level*) 0)
(setf (aref *tags* *level*) 4))
(setf (aref *word* *level*)
(mat* (aref *word* (- *level* 1)) (aref *gens* (aref *tags* *level*))))
(go-forward))))
(compile 'go-forward)
(defun go-backward ()
(decf *level*)
(if (and (/= *level* 0) (not (available-turn-p)))
(go-backward)))
(compile 'go-backward)
(defun available-turn-p ()
(let ((next-next-tag (abs (mod (+ (aref *tags* *level*) 2) 4)))
(pre-tag (- (aref *tags* (+ *level* 1)) 1)))
(if (= next-next-tag 0) (setf next-next-tag 4))
(if (= pre-tag 0) (setf pre-tag 4))
(if (= next-next-tag pre-tag)
nil
t)))
(defun turn-and-go-forward ()
(setf (aref *tags* (+ *level* 1))
(abs (mod (- (aref *tags* (+ *level* 1)) 1) 4)))
(if (= (aref *tags* (+ *level* 1)) 0)
(setf (aref *tags* (+ *level* 1)) 4))
(if (= *level* 0)
(setf (aref *word* 1) (aref *gens* (aref *tags* 1)))
(setf (aref *word* (+ *level* 1))
(mat* (aref *word* *level*)
(aref *gens* (aref *tags* (+ *level* 1))))))
(incf *level*))
(defun branch-termination-p ()
(let ((z (make-array '(4))))
(loop for i from 1 to 3 do
(setf (aref z i)
(mobius-on-point (aref *word* *level*)
(aref *fix-point* (aref *tags* *level*) i))))
(if (or (= *level* *max-level*)
(and (<= (abs (- (aref z 2) (aref z 1))) *epsilon*)
(<= (abs (- (aref z 3) (aref z 2))) *epsilon*)))
(setf *line-points-list*
(append *line-points-list* (list (aref z 1) (aref z 2) (aref z 2) (aref z 3))))
nil)))
(defun explore-limiting-set-with-dfs ()
(labels ((explore ()
(go-forward)
(go-backward)
(turn-and-go-forward)
(if (or (/= *level* 1) (/= (aref *tags* 1) 1))
(explore))))
(explore)))
(compile 'explore-limiting-set-with-dfs)
(defun get-limiting-set (t_a t_b max-level epsilon)
(init max-level epsilon)
(set-gens t_a t_b)
(explore-limiting-set-with-dfs))
(ql:quickload "cl-glut")
(defparameter *width* 500)
(defparameter *height* 500)
(defparameter *magnification* 8)
(defclass my-window (glut:window)
()
(:default-initargs :title "limiting set" :width *width* :height *height*
:mode '(:single :rgb :depth)))
(defmethod glut:display-window :before ((w my-window))
(gl:clear-color 1 1 1 0)
(gl:matrix-mode :projection)
(gl:load-identity)
(gl:ortho 0 *width* *height* 0 -1 1))
(defun set-line-vertexes (lines)
(mapcar #'(lambda (point)
(gl:vertex (* *magnification* (realpart point))
(* *magnification* (imagpart point))
0))
lines))
(defmethod glut:display ((window my-window))
(gl:clear :color-buffer-bit)
(%gl:color-3f 0 0 0)
(gl:push-matrix)
(gl:translate (/ *width* 2) (/ *height* 2) 0)
(gl:begin :lines)
(set-line-vertexes *line-points-list*)
(gl:end)
(gl:pop-matrix)
(gl:flush))
(defun draw-limiting-set ()
(glut:display-window (make-instance 'my-window)))
;; (defun go-forward ()
;; (incf *level*)
;; (setf (aref *tags* *level*)
;; (abs (mod (+ (aref *tags* (- *level* 1)) 1) 4)))
;; (if (= (aref *tags* *level*) 0)
;; (setf (aref *tags* *level*) 4))
;; (setf (aref *word* *level*)
;; (mat* (aref *word* (- *level* 1)) (aref *gens* (aref *tags* *level*)))))
;; (defun go-backward ()
;; (decf *level*))
;; (defun explore-limiting-set-with-dfs ()
;; (labels ((explore ()
;; (loop while (eq (branch-termination-p) nil) do
;; (go-forward))
;; (go-backward)
;; (loop while (and (/= *level* 0) (not (available-turn-p))) do
;; (go-backward))
;; (turn-and-go-forward)))
;; (explore)
;; (loop while (or (/= *level* 1) (/= (aref *tags* 1) 1)) do
;; (explore))))
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
You can’t perform that action at this time.