Library Coq.Numbers.Cyclic.DoubleCyclic.DoubleSqrt
Set Implicit Arguments.
Require Import ZArith.
Require Import BigNumPrelude.
Require Import DoubleType.
Require Import DoubleBase.
Local Open Scope Z_scope.
Section DoubleSqrt.
Variable w :
Type.
Variable w_is_even :
w ->
bool.
Variable w_compare :
w ->
w ->
comparison.
Variable w_0 :
w.
Variable w_1 :
w.
Variable w_Bm1 :
w.
Variable w_WW :
w ->
w ->
zn2z w.
Variable w_W0 :
w ->
zn2z w.
Variable w_0W :
w ->
zn2z w.
Variable w_sub :
w ->
w ->
w.
Variable w_sub_c :
w ->
w ->
carry w.
Variable w_square_c :
w ->
zn2z w.
Variable w_div21 :
w ->
w ->
w ->
w * w.
Variable w_add_mul_div :
w ->
w ->
w ->
w.
Variable w_digits :
positive.
Variable w_zdigits :
w.
Variable ww_zdigits :
zn2z w.
Variable w_add_c :
w ->
w ->
carry w.
Variable w_sqrt2 :
w ->
w ->
w * carry w.
Variable w_pred :
w ->
w.
Variable ww_pred_c :
zn2z w ->
carry (
zn2z w).
Variable ww_pred :
zn2z w ->
zn2z w.
Variable ww_add_c :
zn2z w ->
zn2z w ->
carry (
zn2z w).
Variable ww_add :
zn2z w ->
zn2z w ->
zn2z w.
Variable ww_sub_c :
zn2z w ->
zn2z w ->
carry (
zn2z w).
Variable ww_add_mul_div :
zn2z w ->
zn2z w ->
zn2z w ->
zn2z w.
Variable ww_head0 :
zn2z w ->
zn2z w.
Variable ww_compare :
zn2z w ->
zn2z w ->
comparison.
Variable low :
zn2z w ->
w.
Let wwBm1 :=
ww_Bm1 w_Bm1.
Definition ww_is_even x :=
match x with
|
W0 =>
true
|
WW xh xl =>
w_is_even xl
end.
Let w_div21c x y z :=
match w_compare x z with
|
Eq =>
match w_compare y z with
Eq =>
(C1 w_1, w_0)
|
Gt =>
(C1 w_1, w_sub y z)
|
Lt =>
(C1 w_0, y)
end
|
Gt =>
let x1 :=
w_sub x z in
let (
q,
r) :=
w_div21 x1 y z in
(C1 q, r)
|
Lt =>
let (
q,
r) :=
w_div21 x y z in
(C0 q, r)
end.
Let w_div2s x y s :=
match x with
C1 x1 =>
let x2 :=
w_sub x1 s in
let (
q,
r) :=
w_div21c x2 y s in
match q with
C0 q1 =>
if w_is_even q1 then
(C0 (
w_add_mul_div (
w_pred w_zdigits)
w_1 q1)
, C0 r)
else
(C0 (
w_add_mul_div (
w_pred w_zdigits)
w_1 q1)
, w_add_c r s)
|
C1 q1 =>
if w_is_even q1 then
(C1 (
w_add_mul_div (
w_pred w_zdigits)
w_0 q1)
, C0 r)
else
(C1 (
w_add_mul_div (
w_pred w_zdigits)
w_0 q1)
, w_add_c r s)
end
|
C0 x1 =>
let (
q,
r) :=
w_div21c x1 y s in
match q with
C0 q1 =>
if w_is_even q1 then
(C0 (
w_add_mul_div (
w_pred w_zdigits)
w_0 q1)
, C0 r)
else
(C0 (
w_add_mul_div (
w_pred w_zdigits)
w_0 q1)
, w_add_c r s)
|
C1 q1 =>
if w_is_even q1 then
(C0 (
w_add_mul_div (
w_pred w_zdigits)
w_1 q1)
, C0 r)
else
(C0 (
w_add_mul_div (
w_pred w_zdigits)
w_1 q1)
, w_add_c r s)
end
end.
Definition split x :=
match x with
|
W0 =>
(w_0,w_0)
|
WW h l =>
(h,l)
end.
Definition ww_sqrt2 x y :=
let (
x1,
x2) :=
split x in
let (
y1,
y2) :=
split y in
let (
q,
r) :=
w_sqrt2 x1 x2 in
let (
q1,
r1) :=
w_div2s r y1 q in
match q1 with
C0 q1 =>
let q2 :=
w_square_c q1 in
let a :=
WW q q1 in
match r1 with
C1 r2 =>
match ww_sub_c (
WW r2 y2)
q2 with
C0 r3 =>
(a, C1 r3)
|
C1 r3 =>
(a, C0 r3)
end
|
C0 r2 =>
match ww_sub_c (
WW r2 y2)
q2 with
C0 r3 =>
(a, C0 r3)
|
C1 r3 =>
let a2 :=
ww_add_mul_div (
w_0W w_1)
a W0 in
match ww_pred_c a2 with
C0 a3 =>
(ww_pred a, ww_add_c a3 r3)
|
C1 a3 =>
(ww_pred a, C0 (
ww_add a3 r3)
)
end
end
end
|
C1 q1 =>
let a1 :=
WW q w_Bm1 in
let a2 :=
ww_add_mul_div (
w_0W w_1)
a1 wwBm1 in
(a1, ww_add_c a2 y)
end.
Definition ww_is_zero x :=
match ww_compare W0 x with
Eq =>
true
|
_ =>
false
end.
Definition ww_head1 x :=
let p :=
ww_head0 x in
if (
ww_is_even p)
then p else ww_pred p.
Definition ww_sqrt x :=
if (
ww_is_zero x)
then W0
else
let p :=
ww_head1 x in
match ww_compare p W0 with
|
Gt =>
match ww_add_mul_div p x W0 with
W0 =>
W0
|
WW x1 x2 =>
let (
r,
_) :=
w_sqrt2 x1 x2 in
WW w_0 (
w_add_mul_div
(
w_sub w_zdigits
(
low (
ww_add_mul_div (
ww_pred ww_zdigits)
W0 p)))
w_0 r)
end
|
_ =>
match x with
W0 =>
W0
|
WW x1 x2 =>
WW w_0 (
fst (
w_sqrt2 x1 x2))
end
end.
Variable w_to_Z :
w ->
Z.
Notation wB := (
base w_digits).
Notation wwB := (
base (
ww_digits w_digits)).
Notation "[| x |]" := (
w_to_Z x) (
at level 0,
x at level 99).
Notation "[+| c |]" :=
(
interp_carry 1
wB w_to_Z c) (
at level 0,
x at level 99).
Notation "[-| c |]" :=
(
interp_carry (-1)
wB w_to_Z c) (
at level 0,
x at level 99).
Notation "[[ x ]]" := (
ww_to_Z w_digits w_to_Z x)(
at level 0,
x at level 99).
Notation "[+[ c ]]" :=
(
interp_carry 1
wwB (
ww_to_Z w_digits w_to_Z)
c)
(
at level 0,
x at level 99).
Notation "[-[ c ]]" :=
(
interp_carry (-1)
wwB (
ww_to_Z w_digits w_to_Z)
c)
(
at level 0,
x at level 99).
Notation "[|| x ||]" :=
(
zn2z_to_Z wwB (
ww_to_Z w_digits w_to_Z)
x) (
at level 0,
x at level 99).
Notation "[! n | x !]" := (
double_to_Z w_digits w_to_Z n x)
(
at level 0,
x at level 99).
Variable spec_w_0 :
[|w_0|] = 0.
Variable spec_w_1 :
[|w_1|] = 1.
Variable spec_w_Bm1 :
[|w_Bm1|] = wB - 1.
Variable spec_w_zdigits :
[|w_zdigits|] = Zpos w_digits.
Variable spec_more_than_1_digit: 1
< Zpos w_digits.
Variable spec_ww_zdigits :
[[ww_zdigits]] = Zpos (
xO w_digits).
Variable spec_to_Z :
forall x, 0
<= [|x|] < wB.
Variable spec_to_w_Z :
forall x, 0
<= [[x]] < wwB.
Variable spec_w_WW :
forall h l,
[[w_WW h l]] = [|h|] * wB + [|l|].
Variable spec_w_W0 :
forall h,
[[w_W0 h]] = [|h|] * wB.
Variable spec_w_0W :
forall l,
[[w_0W l]] = [|l|].
Variable spec_w_is_even :
forall x,
if w_is_even x then [|x|] mod 2
= 0
else [|x|] mod 2
= 1.
Variable spec_w_compare :
forall x y,
w_compare x y = Z.compare [|x|] [|y|].
Variable spec_w_sub :
forall x y,
[|w_sub x y|] = ([|x|] - [|y|]) mod wB.
Variable spec_w_square_c :
forall x,
[[ w_square_c x]] = [|x|] * [|x|].
Variable spec_w_div21 :
forall a1 a2 b,
wB/2
<= [|b|] ->
[|a1|] < [|b|] ->
let (
q,
r) :=
w_div21 a1 a2 b in
[|a1|] *wB+ [|a2|] = [|q|] * [|b|] + [|r|] /\
0
<= [|r|] < [|b|].
Variable spec_w_add_mul_div :
forall x y p,
[|p|] <= Zpos w_digits ->
[| w_add_mul_div p x y |] =
([|x|] * (2
^ [|p|]) +
[|y|] / (Z.pow 2 (
(Zpos w_digits) - [|p|])
)) mod wB.
Variable spec_ww_add_mul_div :
forall x y p,
[[p]] <= Zpos (
xO w_digits) ->
[[ ww_add_mul_div p x y ]] =
([[x]] * (2
^ [[p]]) +
[[y]] / (2
^ (Zpos (
xO w_digits)
- [[p]]))) mod wwB.
Variable spec_w_add_c :
forall x y,
[+|w_add_c x y|] = [|x|] + [|y|].
Variable spec_ww_add :
forall x y,
[[ww_add x y]] = ([[x]] + [[y]]) mod wwB.
Variable spec_w_sqrt2 :
forall x y,
wB/ 4
<= [|x|] ->
let (
s,
r) :=
w_sqrt2 x y in
[[WW x y]] = [|s|] ^ 2
+ [+|r|] /\
[+|r|] <= 2
* [|s|].
Variable spec_ww_sub_c :
forall x y,
[-[ww_sub_c x y]] = [[x]] - [[y]].
Variable spec_ww_pred_c :
forall x,
[-[ww_pred_c x]] = [[x]] - 1.
Variable spec_pred :
forall x,
[|w_pred x|] = ([|x|] - 1
) mod wB.
Variable spec_ww_pred :
forall x,
[[ww_pred x]] = ([[x]] - 1
) mod wwB.
Variable spec_ww_add_c :
forall x y,
[+[ww_add_c x y]] = [[x]] + [[y]].
Variable spec_ww_compare :
forall x y,
ww_compare x y = Z.compare [[x]] [[y]].
Variable spec_ww_head0 :
forall x, 0
< [[x]] ->
wwB/ 2
<= 2
^ [[ww_head0 x]] * [[x]] < wwB.
Variable spec_low:
forall x,
[|low x|] = [[x]] mod wB.
Let spec_ww_Bm1 :
[[wwBm1]] = wwB - 1.
Hint Rewrite spec_w_0 spec_w_1 spec_w_WW spec_w_sub
spec_w_add_mul_div spec_ww_Bm1 spec_w_add_c :
w_rewrite.
Lemma spec_ww_is_even :
forall x,
if ww_is_even x then [[x]] mod 2
= 0
else [[x]] mod 2
= 1.
Theorem spec_w_div21c :
forall a1 a2 b,
wB/2
<= [|b|] ->
let (
q,
r) :=
w_div21c a1 a2 b in
[|a1|] * wB + [|a2|] = [+|q|] * [|b|] + [|r|] /\ 0
<= [|r|] < [|b|].
Theorem C0_id:
forall p,
[+|C0 p|] = [|p|].
Theorem add_mult_div_2:
forall w,
[|w_add_mul_div (
w_pred w_zdigits)
w_0 w|] = [|w|] / 2.
Theorem add_mult_div_2_plus_1:
forall w,
[|w_add_mul_div (
w_pred w_zdigits)
w_1 w|] =
[|w|] / 2
+ 2
^ Zpos (
w_digits - 1).
Theorem add_mult_mult_2:
forall w,
[|w_add_mul_div w_1 w w_0|] = 2
* [|w|] mod wB.
Theorem ww_add_mult_mult_2:
forall w,
[[ww_add_mul_div (
w_0W w_1)
w W0]] = 2
* [[w]] mod wwB.
Theorem ww_add_mult_mult_2_plus_1:
forall w,
[[ww_add_mul_div (
w_0W w_1)
w wwBm1]] =
(2
* [[w]] + 1
) mod wwB.
Theorem Zplus_mod_one:
forall a1 b1, 0
< b1 ->
(a1 + b1) mod b1 = a1 mod b1.
Lemma C1_plus_wB:
forall x,
[+|C1 x|] = wB + [|x|].
Theorem spec_w_div2s :
forall a1 a2 b,
wB/2
<= [|b|] ->
[+|a1|] <= 2
* [|b|] ->
let (
q,
r) :=
w_div2s a1 a2 b in
[+|a1|] * wB + [|a2|] = [+|q|] * (2
* [|b|]) + [+|r|] /\ 0
<= [+|r|] < 2
* [|b|].
Theorem wB_div_4: 4
* (wB / 4
) = wB.
Theorem Zsquare_mult:
forall p,
p ^ 2
= p * p.
Theorem Zsquare_pos:
forall p, 0
<= p ^ 2.
Lemma spec_split:
forall x,
[|fst (
split x)
|] * wB + [|snd (
split x)
|] = [[x]].
Theorem mult_wwB:
forall x y,
[|x|] * [|y|] < wwB.
Hint Resolve mult_wwB.
Lemma spec_ww_sqrt2 :
forall x y,
wwB/ 4
<= [[x]] ->
let (
s,
r) :=
ww_sqrt2 x y in
[||WW x y||] = [[s]] ^ 2
+ [+[r]] /\
[+[r]] <= 2
* [[s]].
Lemma spec_ww_is_zero:
forall x,
if ww_is_zero x then [[x]] = 0
else 0
< [[x]].
Lemma wwB_4_2: 2
* (wwB / 4
) = wwB/ 2.
Lemma spec_ww_head1
:
forall x :
zn2z w,
(ww_is_even (
ww_head1 x)
= true) /\
(0
< [[x]] ->
wwB / 4
<= 2
^ [[ww_head1 x]] * [[x]] < wwB).
Theorem wwB_4_wB_4:
wwB / 4
= wB / 4
* wB.
Lemma spec_ww_sqrt :
forall x,
[[ww_sqrt x]] ^ 2
<= [[x]] < ([[ww_sqrt x]] + 1
) ^ 2.
End DoubleSqrt.