proofs/math/typed/f64tools.bend

bend-collections-laws-math@1.0.0.0 · 232 lines · 28 definitions

1import Base
2import bend-collections-laws-math@1.0.0.0/spec/lib/common.bend as C
3import bend-collections-laws-math@1.0.0.0/spec/math/f64.bend as SF
4import bend-collections-laws-math@1.0.0.0/spec/math/w64.bend as SW
5import bend-collections-laws-math@1.0.0.0/src/math/f64.bend as F
6import bend-collections-laws-math@1.0.0.0/src/math/w64.bend as X
7import bend-collections-laws-math@1.0.0.0/src/math/u64.bend as WU
8import bend-collections-laws-math@1.0.0.0/src/math/natural.bend as M
9import bend-collections-laws-math@1.0.0.0/proofs/lib/nat.bend as N
10import bend-collections-laws-math@1.0.0.0/proofs/lib/logic.bend as L
11import bend-collections-laws-math@1.0.0.0/proofs/lib/lemmas/proofs/nat_algebra.bend as NA
12import bend-collections-laws-math@1.0.0.0/proofs/math/typed/width.bend as WW
13import bend-collections-laws-math@1.0.0.0/proofs/math/typed/w64add.bend as WA
14import bend-collections-laws-math@1.0.0.0/proofs/math/typed/w64sh.bend as SH
15import bend-collections-laws-math@1.0.0.0/proofs/math/typed/f64bits.bend as FB
16import bend-collections-laws-math@1.0.0.0/proofs/math/typed/f64light.bend as FL
17import bend-collections-laws-math@1.0.0.0/proofs/math/typed/f64round.bend as FR
18import bend-collections-laws-math@1.0.0.0/proofs/math/typed/f64rtools.bend as RT
19import bend-collections-laws-math@1.0.0.0/proofs/math/typed/f64norm.bend as NM
20import bend-collections-laws-math@1.0.0.0/proofs/math/typed/f64nrp.bend as NR
21import bend-collections-laws-math@1.0.0.0/proofs/math/typed/f64cmp.bend as FC
22import bend-collections-laws-math@1.0.0.0/proofs/math/typed/natcmp.bend as NC
23import bend-collections-laws-math@1.0.0.0/proofs/lib/word.bend as WD
24
25# The two tools of the rounding, conversion and remainder functions: dmant
26# and dexp decode a double into the spec's significand mant (hidden bit
27# included) and scale xexp, and round_w, SoftFloat's normRoundPackToF64 on
28# any 64-bit integer (a sticky shift first when its top bit is set), is the
29# spec's round (Flocq's round_NE of m 2^(x - Z)) for every scale x >= 63.
30
31def v(+x: U32) -> Nat:
32  U32.to_nat(x)
33
34# ---- decoding ----
35
36def dz(+E: Nat, +z: Bool) -> {F.dexp_z(E, z) == SF.pick(Nat, z, Nat.sub(SF.zb(), 1074n), Nat.sub(Nat.add(E, SF.zb()), 1075n)) : Nat}:
37  match z:
38    case True{}:
39      {==}
40    case False{}:
41      Equal.sym(Nat, Nat.sub(Nat.add(E, SF.zb()), 1075n), Nat.add(E, 1925n), Equal.trans(Nat, Nat.sub(Nat.add(E, SF.zb()), 1075n), Nat.add(1925n, E), Nat.add(E, 1925n), FC.esub(E), NA.add_comm(1925n, E)))
42
43def dexp_value(+x: F.F64) -> {F.dexp(x) == SF.xexp(x) : Nat}:
44  match x:
45    case F.Bits{+xl, +xh}:
46      +e1 = Equal.cong(Nat, Nat, t => F.dexp_z(t, Nat.is_eq(t, 0n)), F.exp_field(F.Bits{xl, xh}), SF.efield(F.Bits{xl, xh}), FL.hea(xl, xh))
47      Equal.trans(Nat, F.dexp(F.Bits{xl, xh}), F.dexp_z(SF.efield(F.Bits{xl, xh}), Nat.is_eq(SF.efield(F.Bits{xl, xh}), 0n)), SF.xexp(F.Bits{xl, xh}), e1, dz(SF.efield(F.Bits{xl, xh}), Nat.is_eq(SF.efield(F.Bits{xl, xh}), 0n)))
48
49def dmf(+f: WU.U64) -> {SW.value(F.dmant_z(f, False{})) == SW.value(X.add(f, WU.U64{0, 1048576})) : Nat}:
50  {==}
51
52def dm(+one: Nat, +h1: {one == 1n : Nat}, +f: WU.U64, +Fr: Nat, +hf: {SW.value(f) == Fr : Nat}, +hF: {C.fits(52n, Fr) == True{} : Bool}, +z: Bool) -> {SW.value(F.dmant_z(f, z)) == Nat.add(Fr, C.shift(52n, SF.b2n(Bool.not(z)))) : Nat}:
53  match z:
54    case True{}:
55      +e0 = Equal.trans(Nat, Nat.add(Fr, C.shift(52n, 0n)), Nat.add(Fr, 0n), Fr, Equal.cong(Nat, Nat, t => Nat.add(Fr, t), C.shift(52n, 0n), 0n, WW.shift_zero(52n)), N.add_zero(Fr))
56      Equal.trans(Nat, SW.value(f), Fr, Nat.add(Fr, C.shift(52n, 0n)), hf, Equal.sym(Nat, Nat.add(Fr, C.shift(52n, 0n)), Fr, e0))
57    case False{}:
58      +hFf = L.subst(Nat, t => {C.fits(52n, t) == True{} : Bool}, Fr, SW.value(f), Equal.sym(Nat, SW.value(f), Fr, hf), hF)
59      +e1 = NM.hid(one, h1, f, hFf, 1048576, {==})
60      +e2 = Equal.cong(Nat, Nat, t => Nat.add(t, C.shift(52n, one)), SW.value(f), Fr, hf)
61      # the hidden bit stays `one` / `b2n(not False)`: the checker never
62      # meets the closed 2^52
63      +eo = Equal.trans(Nat, one, 1n, SF.b2n(Bool.not(False{})), h1, {==})
64      +e3 = Equal.cong(Nat, Nat, t => Nat.add(Fr, C.shift(52n, t)), one, SF.b2n(Bool.not(False{})), eo)
65      # the implementation's side unfolds on its own (dmant_z), the spec's
66      # side is the goal's term as written, so the conversion is syntactic
67      +e4 = Equal.trans(Nat, SW.value(X.add(f, WU.U64{0, 1048576})), Nat.add(SW.value(f), C.shift(52n, one)), Nat.add(Fr, C.shift(52n, SF.b2n(Bool.not(False{})))), e1, Equal.trans(Nat, Nat.add(SW.value(f), C.shift(52n, one)), Nat.add(Fr, C.shift(52n, one)), Nat.add(Fr, C.shift(52n, SF.b2n(Bool.not(False{})))), e2, e3))
68      +e5 = dmf(f)
69      Equal.trans(Nat, SW.value(F.dmant_z(f, False{})), SW.value(X.add(f, WU.U64{0, 1048576})), Nat.add(Fr, C.shift(52n, SF.b2n(Bool.not(False{})))), e5, e4)
70
71def dmant_value(+x: F.F64) -> {SW.value(F.dmant(x)) == SF.mant(x) : Nat}:
72  match x:
73    case F.Bits{+xl, +xh}:
74      +e1 = Equal.cong(Nat, Nat, t => SW.value(F.dmant_z(F.frac(F.Bits{xl, xh}), Nat.is_eq(t, 0n))), F.exp_field(F.Bits{xl, xh}), SF.efield(F.Bits{xl, xh}), FL.hea(xl, xh))
75      Equal.trans(Nat, SW.value(F.dmant(F.Bits{xl, xh})), SW.value(F.dmant_z(F.frac(F.Bits{xl, xh}), Nat.is_eq(SF.efield(F.Bits{xl, xh}), 0n))), SF.mant(F.Bits{xl, xh}), e1, dm(1n, {==}, F.frac(F.Bits{xl, xh}), SF.frac(F.Bits{xl, xh}), FL.hfr(xl, xh), FL.hF(xl, xh), Nat.is_eq(SF.efield(F.Bits{xl, xh}), 0n)))
76
77def b2n_fits(+b: Bool) -> {C.fits(1n, SF.b2n(b)) == True{} : Bool}:
78  match b:
79    case True{}:
80      {==}
81    case False{}:
82      {==}
83
84# the significand has at most 53 bits
85def mant_fits(+x: F.F64) -> {C.fits(53n, SF.mant(x)) == True{} : Bool}:
86  match x:
87    case F.Bits{+xl, +xh}:
88      WW.limbs_fit(52n, 1n, SF.frac(F.Bits{xl, xh}), SF.b2n(Bool.not(Nat.is_eq(SF.efield(F.Bits{xl, xh}), 0n))), FL.hF(xl, xh), b2n_fits(Bool.not(Nat.is_eq(SF.efield(F.Bits{xl, xh}), 0n))))
89
90# ---- round_w ----
91
92def zero_eq(+s: Bool) -> {F.zero(s) == SF.zero(s) : F.F64}:
93  match s:
94    case True{}:
95      {==}
96    case False{}:
97      {==}
98
99def nz_mono(+h: Nat, +J: Nat, +hle: {Nat.is_le(h, J) == True{} : Bool}, +hz: {Nat.is_eq(h, 0n) == False{} : Bool}) -> {Nat.is_eq(J, 0n) == False{} : Bool}:
100  match J:
101    case 0n:
102      +e = N.le_antisym(h, 0n, hle, N.zero_le(h))
103      NC.absurd_tf({Nat.is_eq(0n, 0n) == False{} : Bool}, Equal.sym(Bool, True{}, False{}, Equal.trans(Bool, True{}, Nat.is_eq(h, 0n), False{}, Equal.sym(Bool, Nat.is_eq(h, 0n), True{}, L.subst(Nat, t => {Nat.is_eq(t, 0n) == True{} : Bool}, 0n, h, Equal.sym(Nat, h, 0n, e), {==})), hz)))
104    case 1n+j:
105      {==}
106
107# h <= jam(h, l)
108def le_jam(+h: Nat, +l: Nat) -> {Nat.is_le(h, SW.jam(h, l)) == True{} : Bool}:
109  +t = Nat.max(C.bit(h), Nat.min(l, 1n))
110  +ejf = FR.jam_form(h, l)
111  +hhb = WW.le_add_r(C.bit(h), t, Nat.double(C.half(h)), RT.max_ge_l(C.bit(h), Nat.min(l, 1n)))
112  +h1 = L.subst(Nat, z => {Nat.is_le(z, Nat.add(t, Nat.double(C.half(h)))) == True{} : Bool}, Nat.add(C.bit(h), Nat.double(C.half(h))), h, Equal.sym(Nat, h, Nat.add(C.bit(h), Nat.double(C.half(h))), WW.hb(h)), hhb)
113  L.subst(Nat, z => {Nat.is_le(h, z) == True{} : Bool}, Nat.add(t, Nat.double(C.half(h))), SW.jam(h, l), Equal.sym(Nat, SW.jam(h, l), Nat.add(t, Nat.double(C.half(h))), ejf), h1)
114
115# top bit set: a sticky shift by one, then normRoundPack at the next scale
116def rw_hi(+one: Nat, +h1: {one == 1n : Nat}, +s: Bool, +x: Nat, +w: WU.U64, +hx: {Nat.is_le(63n, x) == True{} : Bool}, +hnf: {C.fits(63n, SW.value(w)) == False{} : Bool}) -> {F.rw_top(s, x, w, True{}) == SF.round(s, SW.value(w), x) : F.F64}:
117  +V = SW.value(w)
118  +J = SW.jam(C.high(1n, V), C.low(1n, V))
119  ej1 = SH.shr_jam_value(w, 1n)
120  ej2 = SH.shr_jam_value(w, 1n)
121  ej3 = SH.shr_jam_value(w, 1n)
122  +fh = Equal.trans(Bool, C.fits(63n, C.high(1n, V)), C.fits(64n, V), True{}, Equal.sym(Bool, C.fits(64n, V), C.fits(63n, C.high(1n, V)), FR.fits_hc(1n, 63n, V)), FL.vb64(w))
123  +fJ = Equal.trans(Bool, C.fits(63n, J), C.fits(63n, C.high(1n, V)), True{}, RT.jam_fits(62n, C.high(1n, V), C.low(1n, V)), fh)
124  +nh = FL.nfit_mono(1n, 63n, V, {==}, hnf)
125  +nzJ = nz_mono(C.high(1n, V), J, le_jam(C.high(1n, V), C.low(1n, V)), nh)
126  +h63 = L.subst(Nat, t => {C.fits(63n, t) == True{} : Bool}, J, SW.value(X.shr_jam(w, 1n)), Equal.sym(Nat, SW.value(X.shr_jam(w, 1n)), J, ej1), fJ)
127  +hz = L.subst(Nat, t => {Nat.is_eq(t, 0n) == False{} : Bool}, J, SW.value(X.shr_jam(w, 1n)), Equal.sym(Nat, SW.value(X.shr_jam(w, 1n)), J, ej2), nzJ)
128  +hx1 = N.le_trans(63n, x, Nat.add(x, 1n), hx, N.le_add_right(x, 1n))
129  +ex = NA.add_assoc(x, 1n, 2180n)
130  +r1 = NR.nrp(one, h1, s, Nat.add(x, 2181n), X.shr_jam(w, 1n), Nat.add(x, 1n), ex, hz, h63, hx1)
131  +r2 = Equal.cong(Nat, F.F64, t => SF.round(s, t, Nat.add(x, 1n)), SW.value(X.shr_jam(w, 1n)), J, ej3)
132  +hb = N.le_trans(56n, 64n, M.bit_length(V), {==}, FL.nfit_bl(63n, V, hnf))
133  +r3 = RT.round_jam(s, 1n, V, x, hb)
134  Equal.trans(F.F64, F.norm_round_pack(s, Nat.add(x, 2181n), X.shr_jam(w, 1n)), SF.round(s, SW.value(X.shr_jam(w, 1n)), Nat.add(x, 1n)), SF.round(s, V, x), r1, Equal.trans(F.F64, SF.round(s, SW.value(X.shr_jam(w, 1n)), Nat.add(x, 1n)), SF.round(s, J, Nat.add(x, 1n)), SF.round(s, V, x), r2, Equal.sym(F.F64, SF.round(s, V, x), SF.round(s, J, Nat.add(x, 1n)), r3)))
135
136# the top value 2^63 of the test in rw_z (the word c = 2^31 kept abstract, so
137# the checker never expands 2^31 in unary)
138def top_eq(+one: Nat, +h1: {one == 1n : Nat}, +w: WU.U64, +c: U32, +hc: {c == U32{WD.pw(32n, 31n)} : U32}) -> {X.le(WU.U64{0, c}, w) == Bool.not(C.fits(63n, SW.value(w))) : Bool}:
139  e1 = WA.le_value(WU.U64{0, c}, w)
140  +ec = FB.pwv(31n, {==}, one, h1, c, hc)
141  +e2 = Equal.cong(Nat, Nat, t => Nat.add(0n, C.shift(32n, t)), v(c), C.shift(31n, one), ec)
142  +S = C.shift(32n, C.shift(31n, one))
143  +e30 = Equal.trans(Nat, Nat.add(0n, S), Nat.add(S, 0n), S, NA.add_comm(0n, S), N.add_zero(S))
144  +e3 = Equal.trans(Nat, Nat.add(0n, S), S, C.shift(63n, one), e30, Equal.sym(Nat, C.shift(63n, one), S, WW.shift_comp(32n, 31n, one)))
145  +e4 = Equal.cong(Nat, Bool, t => Nat.is_le(t, SW.value(w)), SW.value(WU.U64{0, c}), C.shift(63n, one), Equal.trans(Nat, SW.value(WU.U64{0, c}), Nat.add(0n, C.shift(32n, C.shift(31n, one))), C.shift(63n, one), e2, e3))
146  Equal.trans(Bool, X.le(WU.U64{0, c}, w), Nat.is_le(SW.value(WU.U64{0, c}), SW.value(w)), Bool.not(C.fits(63n, SW.value(w))), e1, Equal.trans(Bool, Nat.is_le(SW.value(WU.U64{0, c}), SW.value(w)), Nat.is_le(C.shift(63n, one), SW.value(w)), Bool.not(C.fits(63n, SW.value(w))), e4, FB.le_fit(63n, one, h1, SW.value(w))))
147
148def rw_f(+one: Nat, +h1: {one == 1n : Nat}, +s: Bool, +x: Nat, +w: WU.U64, +hx: {Nat.is_le(63n, x) == True{} : Bool}, +hz: {Nat.is_eq(SW.value(w), 0n) == False{} : Bool}, +f: Bool, +hf: {C.fits(63n, SW.value(w)) == f : Bool}) -> {F.rw_top(s, x, w, Bool.not(f)) == SF.round(s, SW.value(w), x) : F.F64}:
149  match f:
150    case True{}:
151      NR.nrp(one, h1, s, Nat.add(x, 2180n), w, x, {==}, hz, hf, hx)
152    case False{}:
153      rw_hi(one, h1, s, x, w, hx, hf)
154
155def rw_c(+one: Nat, +h1: {one == 1n : Nat}, +s: Bool, +x: Nat, +w: WU.U64, +hx: {Nat.is_le(63n, x) == True{} : Bool}, +z: Bool, +hz: {Nat.is_eq(SW.value(w), 0n) == z : Bool}) -> {F.rw_z(s, x, w, z) == SF.round(s, SW.value(w), x) : F.F64}:
156  match z:
157    case True{}:
158      +e1 = Equal.cong(Bool, F.F64, t => SF.pick(F.F64, t, SF.zero(s), SF.round_u(s, SW.value(w), x, Nat.max(Nat.sub(Nat.add(x, M.bit_length(SW.value(w))), 53n), Nat.sub(SF.zb(), 1074n)))), Nat.is_eq(SW.value(w), 0n), True{}, hz)
159      Equal.trans(F.F64, F.zero(s), SF.zero(s), SF.round(s, SW.value(w), x), zero_eq(s), Equal.sym(F.F64, SF.round(s, SW.value(w), x), SF.zero(s), e1))
160    case False{}:
161      +et = Equal.cong(Bool, F.F64, t => F.rw_top(s, x, w, t), X.le(WU.U64{0, 2147483648}, w), Bool.not(C.fits(63n, SW.value(w))), top_eq(one, h1, w, 2147483648, {==}))
162      Equal.trans(F.F64, F.rw_top(s, x, w, X.le(WU.U64{0, 2147483648}, w)), F.rw_top(s, x, w, Bool.not(C.fits(63n, SW.value(w)))), SF.round(s, SW.value(w), x), et, rw_f(one, h1, s, x, w, hx, hz, C.fits(63n, SW.value(w)), {==}))
163
164# round_w is the spec's round
165def rw_value(+s: Bool, +x: Nat, +w: WU.U64, +hx: {Nat.is_le(63n, x) == True{} : Bool}) -> {F.round_w(s, x, w) == SF.round(s, SW.value(w), x) : F.F64}:
166  +ez = Equal.cong(Bool, F.F64, t => F.rw_z(s, x, w, t), X.is_zero(w), Nat.is_eq(SW.value(w), 0n), WA.is_zero_value(w))
167  Equal.trans(F.F64, F.round_w(s, x, w), F.rw_z(s, x, w, Nat.is_eq(SW.value(w), 0n)), SF.round(s, SW.value(w), x), ez, rw_c(1n, {==}, s, x, w, hx, Nat.is_eq(SW.value(w), 0n), {==}))
168
169# ---- the fields as the spec reads them ----
170
171def ea(+x: F.F64) -> {F.exp_field(x) == SF.efield(x) : Nat}:
172  match x:
173    case F.Bits{+xl, +xh}:
174      FL.hea(xl, xh)
175
176def fz(+x: F.F64) -> {X.is_zero(F.frac(x)) == Nat.is_eq(SF.frac(x), 0n) : Bool}:
177  match x:
178    case F.Bits{+xl, +xh}:
179      Equal.trans(Bool, X.is_zero(F.frac(F.Bits{xl, xh})), Nat.is_eq(SW.value(F.frac(F.Bits{xl, xh})), 0n), Nat.is_eq(SF.frac(F.Bits{xl, xh}), 0n), WA.is_zero_value(F.frac(F.Bits{xl, xh})), Equal.cong(Nat, Bool, t => Nat.is_eq(t, 0n), SW.value(F.frac(F.Bits{xl, xh})), SF.frac(F.Bits{xl, xh}), FL.hfr(xl, xh)))
180
181# the scale of a finite double is at least 1926 (the subnormal scale)
182def dge_z(+e: Nat, +z: Bool) -> {Nat.is_le(63n, F.dexp_z(e, z)) == True{} : Bool}:
183  match z:
184    case True{}:
185      {==}
186    case False{}:
187      L.subst(Nat, t => {Nat.is_le(63n, t) == True{} : Bool}, Nat.add(1925n, e), Nat.add(e, 1925n), NA.add_comm(1925n, e), {==})
188
189def dexp_ge(+x: F.F64) -> {Nat.is_le(63n, F.dexp(x)) == True{} : Bool}:
190  dge_z(F.exp_field(x), Nat.is_eq(F.exp_field(x), 0n))
191
192def xexp_ge(+x: F.F64) -> {Nat.is_le(63n, SF.xexp(x)) == True{} : Bool}:
193  L.subst(Nat, t => {Nat.is_le(63n, t) == True{} : Bool}, F.dexp(x), SF.xexp(x), dexp_value(x), dexp_ge(x))
194
195def zero_v(+s: Bool) -> {F.zero(s) == SF.zero(s) : F.F64}:
196  zero_eq(s)
197
198def mw53(+x: F.F64) -> {C.fits(53n, SW.value(F.dmant(x))) == True{} : Bool}:
199  L.subst(Nat, t => {C.fits(53n, t) == True{} : Bool}, SF.mant(x), SW.value(F.dmant(x)), Equal.sym(Nat, SW.value(F.dmant(x)), SF.mant(x), dmant_value(x)), mant_fits(x))
200
201def subsub(+a: Nat, +b: Nat, +h: {Nat.is_le(b, a) == True{} : Bool}) -> {Nat.sub(a, Nat.sub(a, b)) == b : Nat}:
202  +d = Nat.sub(a, b)
203  L.subst(Nat, t => {Nat.sub(t, d) == b : Nat}, Nat.add(b, d), a, N.sub_add(a, b, h), FR.sub_add_l(b, d))
204
205def bzc(+c: Bool) -> {Nat.is_eq(SF.b2n(Bool.not(c)), 0n) == c : Bool}:
206  match c:
207    case True{}:
208      {==}
209    case False{}:
210      {==}
211
212def and_comm(+a: Bool, +b: Bool) -> {Bool.and(a, b) == Bool.and(b, a) : Bool}:
213  match a b:
214    case True{} True{}:
215      {==}
216    case True{} False{}:
217      {==}
218    case False{} True{}:
219      {==}
220    case False{} False{}:
221      {==}
222
223# the significand is zero exactly for the zeros
224def mant_z(+x: F.F64) -> {Nat.is_eq(SF.mant(x), 0n) == SF.is_zero(x) : Bool}:
225  +c = Nat.is_eq(SF.efield(x), 0n)
226  +fz0 = Nat.is_eq(SF.frac(x), 0n)
227  +e1 = FB.z2(SF.frac(x), 52n, SF.b2n(Bool.not(c)))
228  +e2 = Equal.cong(Bool, Bool, b => Bool.and(fz0, b), Nat.is_eq(SF.b2n(Bool.not(c)), 0n), c, bzc(c))
229  Equal.trans(Bool, Nat.is_eq(SF.mant(x), 0n), Bool.and(fz0, Nat.is_eq(SF.b2n(Bool.not(c)), 0n)), SF.is_zero(x), e1, Equal.trans(Bool, Bool.and(fz0, Nat.is_eq(SF.b2n(Bool.not(c)), 0n)), Bool.and(fz0, c), SF.is_zero(x), e2, and_comm(fz0, c)))
230
231def dmz(+x: F.F64) -> {Nat.is_eq(SW.value(F.dmant(x)), 0n) == SF.is_zero(x) : Bool}:
232  Equal.trans(Bool, Nat.is_eq(SW.value(F.dmant(x)), 0n), Nat.is_eq(SF.mant(x), 0n), SF.is_zero(x), Equal.cong(Nat, Bool, t => Nat.is_eq(t, 0n), SW.value(F.dmant(x)), SF.mant(x), dmant_value(x)), mant_z(x))