SRT Square Root Extraction
Function:
(defun quot% (j) (if (zp j) 1 (+ (quot% (1- j)) (* (q% j) (expt (r%) (- j))))))
Function:
(defun rem% (j) (* (expt (r%) j) (- (x%) (* (quot% j) (quot% j)))))
Theorem:
(defthm int-quot-sqrt (implies (natp j) (integerp (* (expt (r%) j) (quot% j)))))
Theorem:
(defthm rem0-sqrt-rewrite (equal (rem% 0) (1- (x%))))
Theorem:
(defthm rem-sqrt-recurrence (implies (natp j) (equal (rem% (+ 1 j)) (- (* (r%) (rem% j)) (* (q% (1+ j)) (+ (* 2 (quot% j)) (* (expt (r%) (- (1+ j))) (q% (1+ j)))))))))
Function:
(defun blo% (j) (+ (* -2 (rho%) (quot% j)) (* (rho%) (rho%) (expt (r%) (- j)))))
Function:
(defun bhi% (j) (+ (* 2 (rho%) (quot% j)) (* (rho%) (rho%) (expt (r%) (- j)))))
Theorem:
(defthm blohi (implies (natp j) (iff (and (<= (expt (- (quot% j) (* (rho%) (expt (r%) (- j)))) 2) (x%)) (>= (expt (+ (quot% j) (* (rho%) (expt (r%) (- j)))) 2) (x%))) (and (<= (blo% j) (rem% j)) (>= (bhi% j) (rem% j))))) :rule-classes nil)
Theorem:
(defthm rem0-sqrt-bounds (and (<= (blo% 0) (rem% 0)) (>= (bhi% 0) (rem% 0))))
Function:
(defun sel-upper-sqrt (k j) (+ (* 2 (+ k (rho%)) (quot% j)) (* (+ k (rho%)) (+ k (rho%)) (expt (r%) (- (1+ j))))))
Function:
(defun sel-lower-sqrt (k j) (+ (* 2 (- k (rho%)) (quot% j)) (* (- k (rho%)) (- k (rho%)) (expt (r%) (- (1+ j))))))
Theorem:
(defthm rem-sqrt-bnds-next (implies (and (natp j) (<= (sel-lower-sqrt (q% (1+ j)) j) (* (r%) (rem% j))) (>= (sel-upper-sqrt (q% (1+ j)) j) (* (r%) (rem% j)))) (and (<= (blo% (1+ j)) (rem% (1+ j))) (>= (bhi% (1+ j)) (rem% (1+ j))))) :rule-classes nil)
Theorem:
(defthm sqrt-containment (implies (natp j) (and (equal (sel-upper-sqrt (a%) j) (* (r%) (bhi% j))) (equal (sel-lower-sqrt (- (a%)) j) (* (r%) (blo% j))))))
Function:
(defun ms4*8 (i j k) (case i (0 (case k (2 12) (1 4) (0 -4) (-1 (if (= j 1) -11 -12)))) (1 (case k (2 (if (= j 2) 15 13)) (1 4) (0 -4) (-1 -13))) (2 (case k (2 15) (1 4) (0 -4) (-1 -15))) (3 (case k (2 16) (1 6) (0 -6) (-1 -16))) (4 (case k (2 18) (1 6) (0 -6) (-1 -18))) (5 (case k (2 20) (1 8) (0 -6) (-1 -20))) (6 (case k (2 20) (1 8) (0 -8) (-1 -20))) (7 (case k (2 22) (1 8) (0 -8) (-1 -22))) (8 (case k (2 24) (1 8) (0 -8) (-1 (if (= j 0) -20 -24))))))
Function:
(defun ms4 (i j k) (/ (ms4*8 i j k) 8))
Function:
(defun i% (j) (* 16 (- (quot% (min (nfix j) 2)) 1/2)))
Function:
(defun select-digit-s4 (a i j) (cond ((<= (ms4 i j 2) a) 2) ((<= (ms4 i j 1) a) 1) ((<= (ms4 i j 0) a) 0) ((<= (ms4 i j -1) a) -1) (t -2)))
Function:
(defun quot%-bnds-inv (j) (and (<= 1/2 (quot% j)) (>= 1 (quot% j))))
Function:
(defun rem%-bnds-inv (j) (and (<= (blo% j) (rem% j)) (>= (bhi% j) (rem% j))))
Function:
(defun approx%-bounds (j k) (and (implies (< (approx% j) (ms4 (i% j) j k)) (< (* 4 (rem% j)) (ms4 (i% j) j k))) (implies (>= (approx% j) (ms4 (i% j) j k)) (> (* 4 (rem% j)) (- (ms4 (i% j) j k) 1/32)))))
Function:
(defun approx%-inv (j) (and (= (q% (1+ j)) (select-digit-s4 (approx% j) (i% j) j)) (approx%-bounds j 2) (approx%-bounds j 1) (approx%-bounds j 0) (approx%-bounds j -1)))
Function:
(defun s4-inv (j) (and (quot%-bnds-inv j) (rem%-bnds-inv j) (approx%-inv j)))
Function:
(defun s4-hyp (j) (if (zp j) (s4-inv 0) (and (s4-inv j) (s4-hyp (1- j)))))
Theorem:
(defthm srt-sqrt-rad-4 (implies (and (= (r%) 4) (= (a%) 2) (natp j) (s4-hyp j)) (and (quot%-bnds-inv (1+ j)) (rem%-bnds-inv (1+ j)))))
Function:
(defun ms8-0 (k) (nth (- 1 k) '(0 -64 -176 -272 -352)))
Function:
(defun ms8-1 (i k) (nth (- 4 k) (nth (/ i 8) '((236 166 96 31 -32 -92 -152 -212) (291 206 121 41 -42 -122 -192 -267) (351 241 141 46 -47 -142 -232 -322) (406 281 171 61 -62 -172 -277 -377) (461 326 191 61 -62 -192 -317 -442)))))
Function:
(defun ms8-2 (i k) (nth (- 4 k) (nth i '((226 161 97 32 -32 -97 -161 -226) (231 165 99 33 -33 -99 -165 -231) (238 170 102 34 -34 -102 -170 -238) (245 175 105 35 -35 -105 -175 -245) (252 180 108 36 -36 -108 -180 -252) (259 185 112 37 -37 -112 -185 -259) (266 190 114 38 -38 -114 -190 -266) (273 195 117 39 -39 -117 -195 -273) (280 200 120 40 -40 -120 -200 -280) (287 205 123 41 -41 -123 -205 -287) (294 210 128 42 -42 -128 -210 -294) (301 215 129 43 -43 -129 -215 -301) (308 220 132 44 -44 -132 -220 -308) (315 225 135 45 -45 -135 -225 -315) (322 230 138 48 -48 -138 -230 -322) (329 235 141 48 -48 -141 -235 -329) (336 240 144 48 -48 -144 -240 -336) (343 245 147 49 -49 -147 -245 -343) (350 250 150 50 -50 -150 -250 -350) (357 255 153 51 -51 -153 -255 -357) (364 260 156 52 -52 -156 -260 -364) (371 265 160 53 -53 -160 -265 -371) (378 270 162 54 -54 -162 -270 -378) (385 275 165 55 -55 -165 -275 -385) (392 280 168 56 -56 -168 -280 -392) (398 285 171 57 -57 -171 -285 -398) (406 290 174 58 -58 -174 -290 -406) (413 295 177 59 -59 -177 -295 -413) (420 300 180 60 -60 -180 -300 -420) (427 305 183 61 -61 -183 -305 -427) (434 310 186 62 -62 -186 -310 -434) (441 315 189 64 -64 -189 -315 -441) (447 319 191 64 -64 -191 -319 -447)))))
Function:
(defun ms8*64 (i j k) (case j (0 (ms8-0 k)) (1 (ms8-1 i k)) (t (ms8-2 i k))))
Function:
(defun ms8 (i j k) (/ (ms8*64 i j k) 64))
Function:
(defun select-digit-s8 (a i j) (cond ((<= (ms8 i j 4) a) 4) ((<= (ms8 i j 3) a) 3) ((<= (ms8 i j 2) a) 2) ((<= (ms8 i j 1) a) 1) ((<= (ms8 i j 0) a) 0) ((<= (ms8 i j -1) a) -1) ((<= (ms8 i j -2) a) -2) ((<= (ms8 i j -3) a) -3) (t -4)))
Function:
(defun i8% (j) (* 64 (- (quot% (min (nfix j) 2)) 1/2)))
Function:
(defun quot%-bnds-inv (j) (and (<= 1/2 (quot% j)) (>= 1 (quot% j))))
Function:
(defun rem%-bnds-inv (j) (and (<= (blo% j) (rem% j)) (>= (bhi% j) (rem% j))))
Function:
(defun approx8%-bounds (j k) (and (implies (< (approx8% j) (ms8 (i8% j) j k)) (< (* 8 (rem% j)) (ms8 (i8% j) j k))) (implies (>= (approx8% j) (ms8 (i8% j) j k)) (> (* 8 (rem% j)) (- (ms8 (i8% j) j k) 1/128)))))
Function:
(defun approx8%-inv (j) (and (= (q% (1+ j)) (select-digit-s8 (approx8% j) (i8% j) j)) (approx8%-bounds j 4) (approx8%-bounds j 3) (approx8%-bounds j 2) (approx8%-bounds j 1) (approx8%-bounds j 0) (approx8%-bounds j -1) (approx8%-bounds j -2) (approx8%-bounds j -3)))
Function:
(defun s8-inv (j) (and (quot%-bnds-inv j) (rem%-bnds-inv j) (approx8%-inv j)))
Function:
(defun s8-hyp (j) (if (zp j) (s8-inv 0) (and (s8-inv j) (s8-hyp (1- j)))))
Theorem:
(defthm srt-sqrt-rad-8 (implies (and (= (r%) 8) (= (a%) 4) (natp j) (s8-hyp j)) (and (quot%-bnds-inv (1+ j)) (rem%-bnds-inv (1+ j)))))