Jest to do zrobienia.

U mnie takie programiki w Auto LISPie robią od ręki studenci, podaję tekst w AutoLispie:
;**************************************************************************
; KZ.LSP 1994
;
; autor: BW - student IV roku
; sprawdzil i nieco zmodyfikowal - MZ
;
;(kolo m z alfan ile) :funkcja rysujaca kolo zebate o srodku w punkcie 0,0
; m - modul kola
; z - liczba zebow kola
; alfan - nominalny kat przyporu
; ile - ilosc odcinkow aproksymujacych ewolwente
;
;(inv x) :inv(x)
; x - kat
;**************************************************************************
(defun inv (x2) (- (/ (sin x2) (cos x2)) x2))
(defun pelna (m rp alfan) (< (- rp m) (* rp (cos alfan))))
(defun *error* (msg) (princ "blad:") (princ msg) (terpri) )
;(trace kolo)
;(trace inv)
;(trace pelna)

(defun kolo (m z alfans ile / rp alfan alfas p1 p2 del1 del2 rprz gam2 fi0 fi1 r0
r1 x0 x1 y0 y1 pkt0 pkt1 delta alfa fipom)
(progn
(setq rp (/ (* m z) 2))
(setq alfan (* alfans (/ pi 180)))
(princ alfans)
(if (pelna m rp alfan)
(setq del1 (inv alfan))
(progn
(setq p1 (/ (- rp m) (* rp (cos alfan))))
(setq alfas (atan (sqrt (- (* p1 p1) 1))))
(setq del1 (- (inv alfan) (inv alfas)))
)
)

(setq p1 (/ (+ rp m) (* rp (cos alfan))))
(setq alfas (atan (sqrt (- (* p1 p1) 1))))
(setq del2 (- (inv alfas) (inv alfan)))
(setq p1 (- rp m) p2 (- rp (* 1.2 m)))
(setq rprz (/ (- (* p1 p1) (* p2 p2)) (* 2.0 p2)))
(setq gam2 (atan (/ rprz p1)))
;
(setq fi0 0.0)
(setq fi1 (- (- (/ pi (* 2.0 z)) del1) gam2))
(setq r0 p2)
(setq r1 p2)
(setq x0 (* r0 (sin fi0)) y0 (* r0 (cos fi0)))
(setq x1 (* r1 (sin fi1)) y1 (* r1 (cos fi1)))
(command "arc" (list x0 y0) "e" (list x1 y1) "r" r0 )
(command "array" "w" (list x1 y1) (list x0 y0) "" "p" (list 0 0) z 360 "y" )
(setq x0 (- 0 x0) x1 (- 0 x1))
(command "arc" (list x1 y1) "e" (list x0 y0) "r" r0 )
(command "array" "w" (list x0 y0) (list x1 y1) "" "p" (list 0 0) z 360 "y" )
;
(setq fi0 fi1)
(setq fi1 (+ fi1 gam2))
(setq r1 (- rp m))
(setq x0 (* r0 (sin fi0)) y0 (* r0 (cos fi0)))
(setq x1 (* r1 (sin fi1)) y1 (* r1 (cos fi1)))
(command "arc" (list x0 y0) "e" (list x1 y1) "r" rprz )
(command "array" "w" (list x1 y1) (list x0 y0) "" "p" (list 0 0) z 360 "y")
(setq x0 (- 0 x0) x1 (- 0 x1))
(command "arc" (list x1 y1) "e" (list x0 y0) "r" rprz )
(command "array" "w" (list x0 y0) (list x1 y1) "" "p" (list 0 0) z 360 "y")

;
(if (pelna m rp alfan)
(progn
(setq fi0 fi1)
(setq r0 (- rp m))
(setq r1 (* rp (cos alfan)))
(setq x0 (* r0 (sin fi0)) y0 (* r0 (cos fi0)))
(setq x1 (* r1 (sin fi1)) y1 (* r1 (cos fi1)))
(command "line" (list x0 y0) (list x1 y1) "")
(command "array" "w" (list x1 y1) (list x0 y0) "" "p" (list 0 0) z 360 "y")
(setq x0 (- 0 x0) x1 (- 0 x1))
(command "line" (list x0 y0) (list x1 y1) "")
(command "array" "w" (list x1 y1) (list x0 y0) "" "p" (list 0 0) z 360 "y")

)
)
;
(setq delta (/ (* m 2.0) ile))
(setq p1 (/ r1 (* rp (cos alfan))))
(setq alfa (atan (sqrt (- (* p1 p1) 1))))
(setq fipom (- fi1 (inv alfa)))
(repeat ile
(progn
(setq r0 r1)
(setq fi0 fi1)
(setq r1 (+ r1 delta))
(setq p1 (/ r1 (* rp (cos alfan))))
(setq alfa (atan (sqrt (- (* p1 p1) 1))))
(setq fi1 (+ fipom (inv alfa)))
(setq x0 (* r0 (sin fi0)) y0 (* r0 (cos fi0)))
(setq x1 (* r1 (sin fi1)) y1 (* r1 (cos fi1)))

(command "line" (list x0 y0) (list x1 y1) "")

(command "array" "w" (list x1 y1) (list x0 y0) "" "p" (list 0 0) z 360 "y")

(setq x0 (- 0 x0) x1 (- 0 x1))

(command "line" (list x0 y0) (list x1 y1) "")

(command "array" "w" (list x1 y1) (list x0 y0) "" "p" (list 0 0) z 360 "y")


)
)
;
(setq r0 r1)
(setq fi0 fi1)
(setq fi1 (/ pi z))
(setq x0 (* r0 (sin fi0)) y0 (* r0 (cos fi0)))
(setq x1 (* r1 (sin fi1)) y1 (* r1 (cos fi1)))

(command "arc" (list x1 y1) "e" (list x0 y0) "r" (+ rp m) )

(command "array" "w" (list x1 y1) (list x0 y0) "" "p" (list 0 0) z 360 "y")

(setq x0 (- 0 x0) x1 (- 0 x1))

(command "arc" (list x0 y0) "e" (list x1 y1) "r" (+ rp m) )

(command "array" "w" (list x0 y0) (list x1 y1) "" "p" (list 0 0) z 360 "y")

)
)
(setq m1 (getreal "\nPodaj modul kola zebatego:"))
(setq z1 (getint "Podaj liczbe zebow: "))
(setq a1 (getreal "Podaj kat przyporu: "))
(setq i1 (getint "Dokladnosc rysowania - podaj liczbe elementow ewolwenty: "))
(setq lim (* 1.5 (/ (* m1 z1) 2.0)))
(setq nlim (- 0 lim))
(setq ldx (list nlim nlim) pg (list lim lim))
(command "limits" ldx pg)
(command "zoom" "all")
(redraw)
(terpri)
(kolo m1 z1 a1 i1)