;;;Pi approxiomation by Archimedes ,then Viete ;;; by Takaya Iwamoto Jan 26,2007 ;;;Ref. "A history of Pi" by Petr Beckmann page 62. ;;; ;;;******************************************************************** ;;; output text file name is pi_number.txt ;;;******************************************************************** ;;; ;;; ;;; (defun c:List_pi_number() (setup_pi_approx) (graph_note) (make_jpg) (setq fil (open "pi_number.txt" "w")) (setq k_times (getint "\nNumber of cycles ?:(def = 17)")) (if (= k_times nil) (setq k_times 17)) (princ (setq title1 "\t\tPi Number List\n\n") fil) ;(princ (setq title2 (strcat "\t " "k" "\t" "lower_pi" "\t" ; "upper_pi" "\n")) fil) (setq n 6 theta (/ pi 6.) nk 0 ) (repeat k_times (setq two_k (* 2 nk) exp_2k (expt 2 nk) theta_2k (/ theta exp_2k) lower_pi (* exp_2k n (sin theta_2k)) upper_pi (* exp_2k n (tan theta_2k)) delta_low (- pi lower_pi) delta_up (- upper_pi pi) y_low (- (+ 12 (log10 delta_low))) y_up (+ 12 (log10 delta_up)) pnt_low (list nk y_low) pnt_up (list nk y_up) app_low (rtos lower_pi 2 16) app_up (rtos upper_pi 2 16) ) (princ "*** ")(princ nk)(princ " : ")(princ delta_low)(princ " : ")(princ delta_up)(terpri) ;(princ y_low)(princ " : ")(princ y_up)(terpri) (princ app_low)(princ " : ")(princ app_up)(terpri) (princ " ")(terpri) ;;Archimedes result (if (= nk 4) (progn (textdisplay (rtos lower_pi 2 6) pnt_low 0.5 0.) (textdisplay (rtos upper_pi 2 6) pnt_up 0.5 0.) ) ) ;;Vieta result (if (and (> nk 7) (= (rem nk 4) 0)) (progn (textdisplay (rtos lower_pi 2 11) pnt_low 0.5 0.) (textdisplay (rtos upper_pi 2 11) pnt_up 0.5 0.) ) ) (if (= nk 0) (progn (setq pnt_low_old pnt_low pnt_up_old pnt_up) (make_pt "0" 1 pnt_low) (make_pt "0" 2 pnt_up) ) (progn (setq pnt_low_new pnt_low pnt_up_new pnt_up) (make_pt "0" 1 pnt_low) (make_pt "0" 2 pnt_up) (make_line_1 "0" 1 pnt_low_old pnt_low_new) (make_line_1 "0" 2 pnt_up_old pnt_up_new) (setq pnt_low_old pnt_low_new pnt_up_old pnt_up_new ) );progn );if ;(princ (setq final_str (strcat "\t" ; nk " , " lower_pi " , " ; upper_pi "\n")) fil) (setq nk (1+ nk)) (make_jpg) (command "_.delay" 1000) );;;repeat loop (make_jpg) ;(final_result) (close fil) );c:List_pi_number ;;; ;;;log with base 10--AutoLISP's log is base "e". ;;; (defun log10( x / inv_log10) (setq inv_log10 (/ 1. (log 10)) ) (princ (* (log x) inv_log10)) );; ;;; ;;; (defun setup_pi_approx() (setup_sysvar) (setvar "PDMODE" 32) (setvar "PDSIZE" -2) (setq pnt_org '(0 0) x_right '(18 0) yaxis_up '(0 12) yaxis_down '(0 -12) upper_right '(20.5 12.5) lower_left '(-3.0 -12.5) ) (make_line_1 "0" 8 pnt_org x_right) (setq horiz_line (entlast)) (make_line_1 "0" 8 yaxis_down yaxis_up) (setq vert_line (entlast)) (command "_.array" vert_line "" "_R" 1 19 1.0) (command "_.array" horiz_line "" "_R" 13 1 1.0) (command "_.array" horiz_line "" "_R" 13 1 -1.0) (make_pt "0" 0 pnt_org) (make_line_1 "0" 0 pnt_org x_right) (make_line_1 "0" 0 yaxis_down yaxis_up) (graph_note) (command "_.zoom" "_W" lower_left upper_right) (command "_.regen") ) ;;; ;;; (defun graph_note() (setq chr_size 0.5 icnt 1) (repeat 12 (setq num (itoa (+ icnt -12)) y_plus (list -1 icnt) y_minus (list -1 (- icnt)) ) (textdisplay num y_plus chr_size 0.) (textdisplay num y_minus chr_size 0.) (setq icnt (1+ icnt)) ) (textdisplay "-12" '(-1 0) chr_size 0.) (textdisplay "Deviation from exact ( 10E-n )" '(-2.5 -9) 1.0 90.) (textdisplay "Lower" '( 5.0 -4) 0.75 0.) (textdisplay "Bound" '( 5.0 -5) 0.75 0.) (textdisplay "Upper" '( 5.0 5.0 ) 0.75 0.) (textdisplay "Bound" '( 5.0 4.0 ) 0.75 0.) (setq icnt 2) (repeat 8 (setq num (itoa icnt) x_loc (list icnt -0.6) ) (textdisplay num x_loc chr_size 0.) (setq icnt (+ 2 icnt)) ) ;;color id (make_Line_1 "0" 3 '(4 -9) '(4 9)) (make_line_1 "0" 4 '(16 -4) '(16 4)) (textdisplay "Archimedes" '(5.5 10.5) 1.0 0.) (textdisplay "Viete" '(14.5 5.0) 1.0 0.) ) ;;; ;;; (defun C:Archimedes_pi() (setup_pi_approx) (setq fil (open "pi_number.txt" "w")) (setq circle_center '(-9 0) circle_radius 5.0 graph_center '(9 0) nk 0 ) (textdisplay "N = " '(-10.3 0.5) 0.7 0.) (textdisplay "< pi <" '(-10.1 -1.0) 0.5 0.) (make_circle_1 "0" 8 circle_center circle_radius) (command "_.vports" 2 "_V") (setvar "CVPORT" 2) (command "_.zoom" "_W" '(-15 -5.5) '(-3 5.5)) (command "_.regen") (setvar "CVPORT" 3) (command "_.zoom" "_W" lower_left upper_right) (setq chr_size 0.5 icnt 1) (graph_note) (make_jpg) ;;;********************************************************************* (while (= (Yes_or_No "Go to next step ?") "_Y") (setq npoly (getint "\nNumber of sides in regular polygon(def =6)")) (if (= npoly nil) (setq npoly 6)) (if (/= nk 0) (progn (entdel in_scribed) (entdel circum_scribed) (entdel npoly_text) (entdel lower_pi_text) (entdel upper_pi_text) ) );;endif (textdisplay (itoa npoly) '(-8.366 0.486) 0.75 0.) (setq npoly_text (entlast)) (setvar "CECOLOR" "1") (command "_.polygon" npoly circle_center "_I" circle_radius) (setq in_scribed (entlast)) (setvar "CECOLOR" "2") (command "_.polygon" npoly circle_center "_C" circle_radius) (setq circum_scribed (entlast)) (command "_.list" in_scribed "") (setq lower_pi (* 0.1 (getvar "PERIMETER"))) (command "_.list" circum_scribed "") (setq upper_pi (* 0.1 (getvar "PERIMETER"))) (setvar "CECOLOR" "0") (textdisplay (rtos lower_pi 2 5) '(-12.77 -1.00) 0.5 0.) (setq lower_pi_text (entlast)) (textdisplay (rtos upper_pi 2 5) '(-8.05 -1.00) 0.5 0.) (setq upper_pi_text (entlast)) (setq result (strcat (itoa npoly) "," (rtos lower_pi 2 8) "," (rtos upper_pi 2 8) )) (princ result)(terpri) (setq delta_low (- pi lower_pi) delta_up (- upper_pi pi) y_low (- (+ 12 (log10 delta_low))) y_up (+ 12 (log10 delta_up)) pnt_low (list nk y_low) pnt_up (list nk y_up) app_low (rtos lower_pi 2 16) app_up (rtos upper_pi 2 16) ) (textdisplay (rtos lower_pi 2 6) pnt_low 0.5 0.) (textdisplay (rtos upper_pi 2 6) pnt_up 0.5 0.) (if (= nk 0) (progn (setq pnt_low_old pnt_low pnt_up_old pnt_up) (make_pt "0" 1 pnt_low) (make_pt "0" 2 pnt_up) ) (progn (setq pnt_low_new pnt_low pnt_up_new pnt_up) (make_pt "0" 1 pnt_low) (make_pt "0" 2 pnt_up) (make_line_1 "0" 1 pnt_low_old pnt_low_new) (make_line_1 "0" 2 pnt_up_old pnt_up_new) (setq pnt_low_old pnt_low_new pnt_up_old pnt_up_new ) );progn );if (setq nk (1+ nk)) (make_jpg) );;while loop (reset_sysvar) ) ;;; ;;; ;;;Huygens ;;; (defun C:Huygens() (setq n_repeat 10 pi_by_6 (/ pi 6.) k_cnt 0 ) (repeat n_repeat (setq c1 (expt 2 k_cnt) c2 (expt 2 (1+ k_cnt)) f1 (* 6 c1) f2 (* 6 c2) theta_1 (/ pi_by_6 c1) theta_2 (/ pi_by_6 c2) sp1 (* f1 (sin theta_1)) sp2 (* f2 (sin theta_2)) bp1 (* f1 (tan theta_1)) bp2 (* f2 (tan theta_2)) x_low (/ (* 3. f1 (sin theta_1)) (+ 2. (cos theta_1))) ) (setq lower_pi (/ (- (* 4. sp2) sp1) 3.) upper_pi (/ (+ (* 2 sp1) bp1) 3.) low_pi (rtos lower_pi 2 13) up_pi (rtos upper_pi 2 13) low_x (rtos x_low 2 13) ) (setq k_cnt (1+ k_cnt)) (princ k_cnt)(terpri) (princ low_pi)(princ " :: ")(princ up_pi)(princ " :: ")(princ low_x)(terpri) );;end of loop ) ;;; ;;; ;;;Snell_Huygens_LB ;;; (defun c:Snell_Huygens_LB( ) (setup_snell_huyg_lb) (reset_sysvar) );; ;;; ;;; (defun setup_snell_huyg_lb() (setup_sysvar) (setvar "PDMODE" 32) (setvar "PDSIZE" -4) ;;lower bound (setq pnt_e '(-2 0) pnt_b '(1 0) pnt_a '(-1 0) pnt_org '(0 0) pnt_f (polar pnt_org (/ pi 3.) 1.0) pnt_tan '(1 0.75) pnt_g (inters pnt_e pnt_f pnt_b pnt_tan nil) dist_bg (distance pnt_b pnt_g) lower_left '(-2.1 -0.5) upper_right '(1.25 1.25) ) ;; (make_pt "0" 0 pnt_org) (make_pt "0" 0 pnt_g) (make_pt "0" 0 pnt_e) (make_pt "0" 0 pnt_a) (make_pt "0" 0 pnt_b) (make_pt "0" 0 pnt_f) (mark_id pnt_org "O" 4 0.1) (mark_id pnt_b "B" 4 0.1) (mark_id pnt_a "A" 4 0.1) (mark_id pnt_e "E" 4 0.1) (mark_id pnt_g "G1" 1 0.1) (mark_id pnt_f "F" 2 0.1) ;; (command "_.arc" "_C" pnt_org pnt_b pnt_e) (make_line_1 "0" 8 pnt_e pnt_b) (make_line_1 "0" 8 pnt_e pnt_g) (make_line_1 "0" 2 pnt_b pnt_g) (make_line_1 "0" 8 pnt_org pnt_f) ;; (textdisplay "r" '(-1.5 0) 0.1 0.) (textdisplay "r" '(-0.5 0) 0.1 0.) (textdisplay "r" '( 0.5 0) 0.1 0.) (textdisplay "arc BF > BG1" '(-0.75 -0.48) 0.2 0.) (textdisplay "Lower" '(-1.68 -0.33) 0.15 0.) (textdisplay "Bound" '(-1.68 -0.57) 0.15 0.) ;; (mark_angle pnt_g pnt_e pnt_b 0.75 "b" 0.15 1) (mark_angle pnt_f pnt_org pnt_b 0.35 "q" 0.15 1) ;; (command "_.zoom" "_W" lower_left upper_right) (command "_.regen") ) ;;;setup_snell_huyg_lb ;;; ;;; ;;;Snell_Huygens_UB ;;; (defun c:Snell_Huygens_UB( ) (setup_snell_huyg_ub) (reset_sysvar) );; ;;; ;;; (defun setup_snell_huyg_ub() (setup_sysvar) (setvar "PDMODE" 32) (setvar "PDSIZE" -4) ;;lower bound (setq pnt_b '(1 0) pnt_a '(-1 0) pnt_org '(0 0) pnt_f (polar pnt_org (/ pi 3.) 1.0) pnt_d (polar pnt_org (/ (* 8. pi) 9.) 1.0) pnt_tan '(1 0.75) pnt_e (inters pnt_d pnt_f pnt_a pnt_b nil) pnt_g (inters pnt_d pnt_f pnt_b pnt_tan nil) dist_bg (distance pnt_b pnt_g) lower_left '(-2.1 -0.5) upper_right '(1.25 1.25) ) ;; (make_pt "0" 0 pnt_org) (make_pt "0" 0 pnt_g) (make_pt "0" 0 pnt_e) (make_pt "0" 0 pnt_a) (make_pt "0" 0 pnt_b) (make_pt "0" 0 pnt_f) (make_pt "0" 0 pnt_d) (mark_id pnt_org "O" 4 0.1) (mark_id pnt_b "B" 4 0.1) (mark_id pnt_a "A" 4 0.1) (mark_id pnt_e "E" 4 0.1) (mark_id pnt_g "G2" 1 0.1) (mark_id pnt_f "F" 2 0.1) (mark_id pnt_d "D" 3 0.1) ;; (command "_.arc" "_C" pnt_org pnt_b pnt_e) (make_line_1 "0" 8 pnt_e pnt_b) (make_line_1 "0" 8 pnt_e pnt_g) (make_line_1 "0" 2 pnt_b pnt_g) (make_line_1 "0" 8 pnt_org pnt_f) (make_line_1 "0" 8 pnt_org pnt_d) ;; (textdisplay "r" '(-1.5 0.21) 0.1 0.) (textdisplay "r" '(-0.47 0.23) 0.1 0.) (textdisplay "r" '( 0.5 0) 0.1 0.) (textdisplay "arc BF < BG2" '(-0.75 -0.48) 0.2 0.) (textdisplay "Upper" '(-1.68 -0.33) 0.15 0.) (textdisplay "Bound" '(-1.68 -0.57) 0.15 0.) ;; (mark_angle pnt_g pnt_e pnt_b 0.50 "q/3" 0.15 1) (mark_angle pnt_f pnt_org pnt_b 0.35 "q" 0.15 1) ;; (command "_.zoom" "_W" lower_left upper_right) (command "_.regen") ) ;;;setup_snell_huyg_ub ;;; ;;;SNELL_PI ;;; (defun c:Snell_pi() (setq k_times (getint "\nNumber of cycles ?:(def = 10)")) (if (= k_times nil) (setq k_times 10)) (setq hexa (/ pi 6.) nk 0 ) (repeat k_times (setq two_k (* 2 nk) exp_2k (expt 2 nk) theta (/ hexa exp_2k) alpha (/ theta 3.) lbf (* 18. exp_2k) ubf (* 6. exp_2k) lb_num (sin theta) lb_denom (+ 2 (cos theta)) ub_1 (+ (* 2. (cos alpha)) 1) ub_2 (tan alpha) lower_pi (/ (* lbf lb_num) lb_denom) upper_pi (* ubf ub_1 ub_2) delta_low (- pi lower_pi) delta_up (- upper_pi pi) app_low (rtos lower_pi 2 16) app_up (rtos upper_pi 2 16) ) (princ "** ")(princ nk)(princ " **")(terpri) (princ app_low)(princ " : ")(princ app_up)(terpri) (princ " : ")(princ delta_low)(princ " : ")(princ delta_up)(terpri) (princ " ")(terpri) (setq nk (1+ nk)) ) )