#|
Copyright (C) 1994 by Yuan Yu.  All Rights Reserved.
This script is hereby placed in the public domain, and therefore unlimited
editing and redistribution is permitted.
NO WARRANTY
Yuan Yu PROVIDES ABSOLUTELY NO WARRANTY.  THE EVENT SCRIPT IS PROVIDED "AS IS"
WITHOUT WARRANTY OF ANY KIND, EITHER EXPRESS OR IMPLIED, INCLUDING, BUT NOT
LIMITED TO, ANY IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A
PARTICULAR PURPOSE.  THE ENTIRE RISK AS TO THE QUALITY AND PERFORMANCE OF THE
SCRIPT IS WITH YOU.  SHOULD THE SCRIPT PROVE DEFECTIVE, YOU ASSUME THE COST OF
ALL NECESSARY SERVICING, REPAIR OR CORRECTION.
IN NO EVENT WILL Yuan Yu BE LIABLE TO YOU FOR ANY DAMAGES, ANY LOST PROFITS,
LOST MONIES, OR OTHER SPECIAL, INCIDENTAL OR CONSEQUENTIAL DAMAGES ARISING OUT
OF THE USE OR INABILITY TO USE THIS SCRIPT (INCLUDING BUT NOT LIMITED TO LOSS
OF DATA OR DATA BEING RENDERED INACCURATE OR LOSSES SUSTAINED BY THIRD
PARTIES), EVEN IF YOU HAVE ADVISED US OF THE POSSIBILITY OF SUCH DAMAGES, OR
FOR ANY CLAIM BY ANY OTHER PARTY.
|#
;       Proof of the Correctness of an Integer Square Root Program
;

Event: Start with the library "mc20-2" using the compiled version.

#|
This is a revisit to our ISQRT proof.  Dr. Don Good talked with
Dr. Steve Zeigler, vice president for Ada products at Verdix, and
they agreed to let me publicize the machine codes generated by
their present Ada compiler.
The following Ada function ISQRT computes the integer square root of a
given nonnegative integer i.  This is our third proof about ISQRT.
function isqrt (i:integer) return integer is
 j : integer := (i / 2);
begin
 while ((i / j) < j) loop
   j := (j + (i / j)) / 2;

end loop;
return j;
end isqrt;
Here is the MC68020 assembly code of the above ISQRT program.  The code
is from Dr. Steve Zeigler and generated by their present Ada compiler.
1  function isqrt (i:integer) return integer is
00000: link.w      a6, #-04
2     j : integer := (i / 2);
00004: move.l      d2, d1
00006: bge.b       06    -> 0e
00008: addi.l      #01, d1
0000e: asr.l       #01, d1
3  begin
4     while not ((i / j) >= j) loop
00010: move.l      d2, d0
00012: divsl.l     d1, d0:d0
00016: trapv
00018: cmp.l       d0, d1
0001a: ble.b       01c    -> 038
5         j := (j + (i / j)) / 2;
0001c: add.l       d1, d0
0001e: trapv
00020: move.l      d0, d1
00022: bge.b       06    -> 02a
00024: addi.l      #01, d1
0002a: asr.l       #01, d1
4     while not ((i / j) >= j) loop
0002c: move.l      d2, d0
0002e: divsl.l     d1, d0:d0
00032: trapv
6  end loop;
00034: cmp.l       d0, d1
00036: bgt.b       -01c    -> 01c
7  return j;
00038: move.l      d1, d0
0003a: unlk        a6
0003c: rts
8 end isqrt;
0x4e56  0xfffc  0x2202  0x6c06  0x0681  0x0000  0x0001  0xe281
0x2002  0x4c41  0x0800  0x4e76  0xb280  0x6f1c  0xd081  0x4e76
0x2200  0x6c06  0x0681  0x0000  0x0001  0xe281  0x2002  0x4c41

0x0800  0x4e76  0xb280  0x6ee4  0x2001  0x4e5e  0x4e75
'(78      86      255     252     34      2       108     6
  6       129     0       0       0       1       226     129
  32      2       76      65      8       0       78      118
  178     128     111     28      208     129     78      118
  34      0       108     6       6       129     0       0
  0       1       226     129     32      2       76      65
  8       0       78      118     178     128     110     228
  32      1       78      94      78      117)
|#
; in the logic, the above program is defined by (isqrt-code).

Definition: isqrt-code

= '(78 86 255 252 34 2 108 6 6 129 0 0 0 1 226 129 32 2
      76 65 8 0 78 118 178 128 111 28 208 129 78 118 34 0
      108 6 6 129 0 0 0 1 226 129 32 2 76 65 8 0 78 118
      178 128 110 228 32 1 78 94 78 117)

Definition: sq (x) = (x ∗ x)

; isqrt1 is a function in the Logic simulating the loop of the above code.

Definition: isqrt1 (i, j ) = if j ≃ 0 then fix (i)

elseif (i ÷ j ) < j then isqrt1 (i, (j + (i ÷ j )) ÷ 2)

else fix (j ) endif

; isqrt specifies the semantics of ISQRT in the Logic.  To see why
; the function isqrt computes the square root for any nonnegative
; integer input, please refer to the proof in file isqrt.events.

Definition: isqrt (i) = let j1 be ((i ÷ 2) + (i ÷ (i ÷ 2))) ÷ 2

in

if i < sq (i ÷ 2) then isqrt1 (i, j1 )

else i ÷ 2 endif endlet

; the computation time of this program.

Definition: isqrt1-t (i, j ) = if j ≃ 0 then 0

elseif (i ÷ j ) < j

then splus (10, isqrt1-t (i, (j + (i ÷ j )) ÷ 2)) else 8 endif

Definition: isqrt-t (i) = let j1 be ((i ÷ 2) + (i ÷ (i ÷ 2))) ÷ 2

in

if i < sq (i ÷ 2) then splus (14, isqrt1-t (i, j1 )) else 12 endif endlet

; enable a few functions.

Event: Enable iplus.

Event: Enable integerp.

Event: Enable iquotient.

Event: Enable ilessp.

; disable a few functions.

Event: Disable remainder.

Event: Disable quotient.

Theorem: isqrt-no-overflow (int-rangep (2 ∗ j , n) ∧ ((i ÷ j ) < j )) → int-rangep (j + (i ÷ j ), n)

Theorem: j-nonzerop ((1 < i) ∧ (0 < j )) → (((j + (i ÷ j )) ÷ 2) ≠ 0)

Theorem: j-int-rangep (int-rangep (2 ∗ j , n) ∧ ((i ÷ j ) < j )) → int-rangep (2 ∗ ((j + (i ÷ j )) ÷ 2), n)

; an induction hint.

Definition: isqrt-induct (s, i, j ) = if j ≃ 0 then t

elseif (i ÷ j ) < j then isqrt-induct (stepn (s, 10), i, (j + (i ÷ j )) ÷ 2) else t endif

; the properties of the initial state.

Definition: isqrt-statep (s, i) = ((mc-status (s) = 'running)

∧ evenp (mc-pc (s)) ∧ rom-addrp (mc-pc (s), mc-mem (s), 70) ∧ mcode-addrp (mc-pc (s), mc-mem (s), isqrt-code) ∧ ram-addrp (sub (32, 8, read-sp (s)), mc-mem (s), 12) ∧ (i = iread-dn (32, 2, s)) ∧ ilessp (1, i))

; an intermediate state s0.

Definition: isqrt-s0p (s, i, j ) = ((mc-status (s) = 'running)

∧ evenp (mc-pc (s)) ∧ rom-addrp (sub (32, 44, mc-pc (s)), mc-mem (s), 70) ∧ mcode-addrp (sub (32, 44, mc-pc (s)), mc-mem (s), isqrt-code) ∧ ram-addrp (sub (32, 4, read-an (32, 6, s)), mc-mem (s), 12) ∧ (i = iread-dn (32, 2, s)) ∧ (j = iread-dn (32, 1, s)) ∧ int-rangep (2 ∗ j , 32) ∧ ilessp (1, i) ∧ ilessp (0, j ))

Theorem: initial-j-int-rangep int-rangep (i, n) → int-rangep (2 ∗ (i ÷ 2), n)

; from the initial state to exit.

Theorem: isqrt-s-exit (isqrt-statep (s, i) ∧ (i ≮ ((i ÷ 2) ∗ (i ÷ 2)))) → ((mc-status (stepn (s, 12)) = 'running)

∧ (mc-pc (stepn (s, 12)) = rts-addr (s)) ∧ (iread-dn (32, 0, stepn (s, 12)) = (i ÷ 2)) ∧ (read-rn (32, 14, mc-rfile (stepn (s, 12)))

= read-rn (32, 14, mc-rfile (s))) ∧ (read-rn (32, 15, mc-rfile (stepn (s, 12))) = add (32, read-sp (s), 4)))

Theorem: isqrt-s-exit-rfile (isqrt-statep (s, i) ∧ (i ≮ ((i ÷ 2) ∗ (i ÷ 2))) ∧ d2-7a2-5p (rn)) → (read-rn (oplen, rn, mc-rfile (stepn (s, 12)))

= read-rn (oplen, rn, mc-rfile (s)))

Theorem: isqrt-s-exit-mem (isqrt-statep (s, i) ∧ (i ≮ ((i ÷ 2) ∗ (i ÷ 2))) ∧ disjoint (x, k, sub (32, 8, read-sp (s)), 12)) → (read-mem (x, mc-mem (stepn (s, 12)), k) = read-mem (x, mc-mem (s), k))

; from the initial state to s0.

Theorem: isqrt-s-s0 let j1 be ((i ÷ 2) + (i ÷ (i ÷ 2))) ÷ 2 in (isqrt-statep (s, i) ∧ (i < ((i ÷ 2) ∗ (i ÷ 2)))) → (isqrt-s0p (stepn (s, 14), i, j1 )

∧ (linked-rts-addr (stepn (s, 14)) = rts-addr (s)) ∧ (linked-a6 (stepn (s, 14)) = read-an (32, 6, s)) ∧ (read-rn (32, 14, mc-rfile (stepn (s, 14))) = sub (32, 4, read-sp (s)))) endlet

Theorem: isqrt-s-s0-rfile (isqrt-statep (s, i) ∧ (i < ((i ÷ 2) ∗ (i ÷ 2))) ∧ d2-7a2-5p (rn)) → (read-rn (oplen, rn, mc-rfile (stepn (s, 14)))

= read-rn (oplen, rn, mc-rfile (s)))

Theorem: isqrt-s-s0-mem (isqrt-statep (s, i) ∧ (i < ((i ÷ 2) ∗ (i ÷ 2))) ∧ disjoint (x, k, sub (32, 8, read-sp (s)), 12)) → (read-mem (x, mc-mem (stepn (s, 14)), k) = read-mem (x, mc-mem (s), k))

; from s0 to exit (base case), or from s0 to s0 (induction case).
; base case: s0 --> exit.

Theorem: isqrt-s0-exit-base (isqrt-s0p (s, i, j ) ∧ ((i ÷ j ) ≮ j )) → ((mc-status (stepn (s, 8)) = 'running)

∧ (mc-pc (stepn (s, 8)) = linked-rts-addr (s))

∧ (iread-dn (32, 0, stepn (s, 8)) = j )

∧ (read-rn (32, 14, mc-rfile (stepn (s, 8))) = linked-a6 (s)) ∧ (read-rn (32, 15, mc-rfile (stepn (s, 8))) = add (32, read-an (32, 6, s), 8)))

Theorem: isqrt1-s0-exit-rfile-base (isqrt-s0p (s, i, j ) ∧ ((i ÷ j ) ≮ j ) ∧ d2-7a2-5p (rn)) → (read-rn (oplen, rn, mc-rfile (stepn (s, 8)))

= read-rn (oplen, rn, mc-rfile (s)))

Theorem: isqrt1-s0-exit-mem-base (isqrt-s0p (s, i, j ) ∧ ((i ÷ j ) ≮ j )) → (read-mem (x, mc-mem (stepn (s, 8)), k) = read-mem (x, mc-mem (s), k))

; induction case: s0 --> s0.

Theorem: isqrt-s0-s0 (isqrt-s0p (s, i, j ) ∧ ((i ÷ j ) < j )) → (isqrt-s0p (stepn (s, 10), i, (j + (i ÷ j )) ÷ 2)

∧ (read-rn (oplen, 14, mc-rfile (stepn (s, 10))) = read-rn (oplen, 14, mc-rfile (s))) ∧ (linked-a6 (stepn (s, 10)) = linked-a6 (s)) ∧ (linked-rts-addr (stepn (s, 10)) = linked-rts-addr (s)))

Theorem: isqrt-s0-s0-rfile (isqrt-s0p (s, i, j ) ∧ ((i ÷ j ) < j ) ∧ d2-7a2-5p (rn)) → (read-rn (oplen, rn, mc-rfile (stepn (s, 10)))

= read-rn (oplen, rn, mc-rfile (s)))

Theorem: isqrt-s0-s0-mem (isqrt-s0p (s, i, j ) ∧ ((i ÷ j ) < j ) ∧ disjoint (x, k, sub (32, 8, read-an (32, 6, s)), 20)) → (read-mem (x, mc-mem (stepn (s, 10)), k) = read-mem (x, mc-mem (s), k))

Theorem: isqrt-s0p-j nonzerop isqrt-s0p (s, i, j ) → ((j ≠ 0) ∧ (j ∈ N))

Event: Disable isqrt-statep.

Event: Disable isqrt-s0p.

; put together:  s0 --> exit.

Theorem: isqrt-s0-exit let sn be stepn (s, isqrt1-t (i, j )) in isqrt-s0p (s, i, j ) → ((mc-status (sn) = 'running)

∧ (mc-pc (sn) = linked-rts-addr (s)) ∧ (iread-dn (32, 0, sn) = isqrt1 (i, j )) ∧ (read-rn (32, 14, mc-rfile (sn)) = linked-a6 (s)) ∧ (read-rn (32, 15, mc-rfile (sn))

= add (32, read-an (32, 6, s), 8))) endlet

Theorem: isqrt-s0-exit-rfile (isqrt-s0p (s, i, j ) ∧ d2-7a2-5p (rn)) → (read-rn (oplen, rn, mc-rfile (stepn (s, isqrt1-t (i, j ))))

= read-rn (oplen, rn, mc-rfile (s)))

Theorem: isqrt-s0-exit-mem (isqrt-s0p (s, i, j ) ∧ disjoint (x, k, sub (32, 8, read-an (32, 6, s)), 20)) → (read-mem (x, mc-mem (stepn (s, isqrt1-t (i, j ))), k)

= read-mem (x, mc-mem (s), k))

Event: Disable isqrt-s0p-j nonzerop.

; ISQRT program is correct.
; after an execution of this program,  the machine state satisfies:
;  0. normal exit.
;  1. the pc is returned to the next instruction of the caller.
;  2. the result -- ISQRT(i), is stored in D0.
;  3. a6, used by LINK, is restored to its original content.
;  4. the stack pointer sp(a7) is updated correctly to pop off one frame.
; the registers not touched by this program still have their original values.
; the memory is correctly changed -- the local effect on memory.

Theorem: isqrt-correctness let sn be stepn (s, isqrt-t (i)) in isqrt-statep (s, i) → ((mc-status (sn) = 'running)

∧ (mc-pc (sn) = rts-addr (s))

∧ (read-rn (32, 14, mc-rfile (sn)) = read-an (32, 6, s)) ∧ (read-rn (32, 15, mc-rfile (sn))

= add (32, read-an (32, 7, s), 4)) ∧ (d2-7a2-5p (rn)

→ (read-rn (oplen, rn, mc-rfile (sn))

= read-rn (oplen, rn, mc-rfile (s))))

∧ (disjoint (x, k, sub (32, 12, read-sp (s)), 20) → (read-mem (x, mc-mem (sn), k)

= read-mem (x, mc-mem (s), k)))

∧ (iread-dn (32, 0, sn) = isqrt (i))) endlet

Event: Disable isqrt-t.

; start to prove isqrt correct.

Theorem: isqrt1-lower-bound (j ≄ 0) → (i ≮ sq (isqrt1 (i, j )))

Theorem: quotient-by-2 ((x ÷ 2) + (x ÷ 2)) ≮ (x − 1)

Theorem: main-trick sq (1 + ((j + k) ÷ 2)) ≮ ((j ∗ k) + j )

Theorem: sq-add1-non-zero sq (1 + x) ≠ 0

Theorem: main (j ≄ 0) → (i < sq (1 + ((j + (i ÷ j )) ÷ 2)))

Theorem: isqrt1-upper-bound (i < sq (1 + j )) → (i < sq (1 + isqrt1 (i, j )))

Theorem: isqrt->isqrt1 (1 < i) → (i < sq (1 + (i ÷ 2)))

; (isqrt i) is the square root of i:  (isqrt i)^2 <= i < [(isqrt i)+1]^2.

Theorem: isqrt-logic-correctness (1 < i) → ((i < sq (1 + isqrt (i))) ∧ (i ≮ sq (isqrt (i))))

; a simple time analysis.

Theorem: quotient-mono-1 ((y ≤ z) ∧ (y ≄ 0) ∧ (z ≄ 0)) → ((x ÷ y) ≮ (x ÷ z))

Theorem: mean-lessp-1 (x ≮ y) → (x ≮ ((x + y) ÷ 2))

Theorem: isqrt1-t-la-0 let j1 be (j + (i ÷ j )) ÷ 2 in (((i ÷ j ) ≤ j ) ∧ (1 < i) ∧ (0 < j )) → ((x − (i ÷ j )) ≮ (x − (i ÷ j1 ))) endlet

Theorem: mean-difference-1-1 (i ≤ j ) → ((((j + i) ÷ 2) − i) = ((j − i) ÷ 2))

Theorem: isqrt1-t-la-1 let j1 be (j + (i ÷ j )) ÷ 2 in (((i ÷ j ) < j ) ∧ (1 < i) ∧ (0 < j )) → (((j − (i ÷ j )) ÷ 2) ≮ (j1 − (i ÷ j1 ))) endlet

Theorem: isqrt1-t-la-2 let j1 be (j + (i ÷ j )) ÷ 2 in (((i ÷ j ) < j ) ∧ (1 < i) ∧ (0 < j )) → (log (2, (j − (i ÷ j )) ÷ 2) ≮ log (2, j1 − (i ÷ j1 ))) endlet

Theorem: isqrt-t1-ubound ((1 < i) ∧ (0 < j )) → ((18 + (10 ∗ log (2, j − (i ÷ j )))) ≮ isqrt1-t (i, j ))

Theorem: isqrt-t->isqrt1-t (1 < i) → ((4 + isqrt1-t (i, i ÷ 2)) ≮ isqrt-t (i))

Theorem: log-mono log (b, y) ≮ log (b, y − x)

Theorem: isqrt-t-ubound-la (1 < i) → ((12 + (10 ∗ log (2, i))) ≮ isqrt-t (i))

Theorem: isqrt-t-ubound ((i < exp (2, 31)) ∧ (1 < i)) → (isqrt-t (i) ≤ 322)

Index

add, 6–8

d2-7a2-5p, 6–8 disjoint, 6–9

evenp, 5 exp, 10

ilessp, 5 initial-j-int-rangep, 5 int-rangep, 4, 5 iread-dn, 5, 7–9 isqrt, 3, 9 isqrt->isqrt1, 9 isqrt-code, 3, 5 isqrt-correctness, 8 isqrt-induct, 5 isqrt-logic-correctness, 9 isqrt-no-overflow, 4 isqrt-s-exit, 5 isqrt-s-exit-mem, 6 isqrt-s-exit-rfile, 6 isqrt-s-s0, 6 isqrt-s-s0-mem, 6 isqrt-s-s0-rfile, 6 isqrt-s0-exit, 8 isqrt-s0-exit-base, 6 isqrt-s0-exit-mem, 8 isqrt-s0-exit-rfile, 8 isqrt-s0-s0, 7 isqrt-s0-s0-mem, 7 isqrt-s0-s0-rfile, 7 isqrt-s0p, 5–8 isqrt-s0p-j nonzerop, 7 isqrt-statep, 5, 6, 8 isqrt-t, 4, 8, 10 isqrt-t->isqrt1-t, 10 isqrt-t-ubound, 10 isqrt-t-ubound-la, 10 isqrt-t1-ubound, 10 isqrt1, 3, 8, 9 isqrt1-lower-bound, 9 isqrt1-s0-exit-mem-base, 7 isqrt1-s0-exit-rfile-base, 7 isqrt1-t, 3, 4, 8, 10 isqrt1-t-la-0, 10 isqrt1-t-la-1, 10 isqrt1-t-la-2, 10 isqrt1-upper-bound, 9

j-int-rangep, 4 j-nonzerop, 4

linked-a6, 6–8

linked-rts-addr, 6–8

log, 10

log-mono, 10

main, 9 main-trick, 9 mc-mem, 5–9 mc-pc, 5, 6, 8 mc-rfile, 5–9 mc-status, 5, 6, 8 mcode-addrp, 5 mean-difference-1-1, 10 mean-lessp-1, 9

quotient-by-2, 9 quotient-mono-1, 9

ram-addrp, 5 read-an, 5–8 read-mem, 6–9 read-rn, 5–9 read-sp, 5, 6, 9 rom-addrp, 5 rts-addr, 5, 6, 8

splus, 4 sq, 3, 4, 9 sq-add1-non-zero, 9 stepn, 5–8 sub, 5–9