-
-
Save ksk/df4123b361a0faf2e5c41e7729d61f25 to your computer and use it in GitHub Desktop.
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Require Import List PeanoNat Arith Morphisms Setoid Recdef Lia. | |
Import ListNotations. | |
Set Implicit Arguments. | |
Ltac optimize db := | |
eexists; intros; | |
match goal with | |
|- ?X = ?Y => | |
let P := fresh in | |
set (P := fun y => y = Y); | |
enough (P X) by auto | |
end; | |
repeat autounfold with db; | |
autorewrite with db; | |
reflexivity. | |
Definition sqr_list xs := map (fun x => x * x) xs. | |
Definition sqr_add2_list xs := map (fun x => x + 2) (sqr_list xs). | |
Hint Unfold sqr_list sqr_add2_list : my_db. | |
Hint Rewrite map_map : my_db. | |
Definition sa_optimized_sig : { f | forall xs, sqr_add2_list xs = f xs } := | |
ltac:(optimize my_db). | |
Eval simpl in proj1_sig sa_optimized_sig. | |
(* result : | |
[[[ | |
= map (fun x : nat => x * x + 2) | |
: list nat -> list nat | |
]]] | |
*) | |
Lemma fold_left_map : forall (A B C : Type)(f : A -> B) g l (init:C), | |
fold_left g (map f l) init = fold_left (fun x y => g x (f y)) l init. | |
Proof. | |
induction l; simpl; auto. | |
Qed. | |
Definition sum xs := fold_left Nat.add xs 0. | |
Definition sum_of_square xs := sum (map (fun x => x * x) xs). | |
Hint Unfold sum sum_of_square : my_db. | |
Hint Rewrite fold_left_map map_map : my_db. | |
Definition ss_optimized_sig : { prog | forall xs, sum_of_square xs = prog xs } | |
:= ltac:(optimize my_db). | |
Eval simpl in proj1_sig ss_optimized_sig. | |
(* result : | |
[[[ | |
= fun xs : list nat => fold_left (fun x y : nat => x + y * y) xs 0 | |
: list nat -> nat | |
]]] | |
*) | |
Require Import ZArith. | |
Open Scope Z_scope. | |
Fixpoint divide_list (A : Type)(xs : list A) := | |
match xs with | |
| [] => ([],[]) | |
| x::xs' => | |
let (xs0,xs1) := divide_list xs' in | |
(x::xs1,xs0) | |
end. | |
Fixpoint map2 (A B C : Type)(f : A -> B -> C)(xs : list A)(ys : list B) := | |
match xs,ys with | |
| [],_ => [] | |
| _,[] => [] | |
| x::xs',y::ys' => f x y :: map2 f xs' ys' | |
end. | |
Section Fft. | |
Variable N : Z. | |
Variable Npos: 0 <= N. | |
Definition modulus := 2 ^ (2 ^ N) + 1. | |
Definition add x y := | |
(x + y) mod modulus. | |
Definition sub x y := | |
(x - y) mod modulus. | |
Definition mul x y := | |
(x * y) mod modulus. | |
Definition pow_pos x := | |
Pos.iter (mul x) 1. | |
Definition pow x y := | |
match y with | |
| 0 => 1 | |
| Z.pos p => pow_pos x p | |
| Z.neg _ => 0 | |
end. | |
Definition prt n : Z := | |
pow 2 (pow 2 (N + 1 - n)). | |
Definition pow_neg1 (n : Z) := | |
if Z.even n then 1 else -1. | |
Definition shiftl x y := | |
Z.of_nat (Nat.shiftl (Z.to_nat (x mod modulus)) (Z.to_nat y)) mod modulus. | |
Lemma modulus_gt1: 1 < modulus. | |
Proof. | |
unfold modulus. | |
now apply Z.lt_add_pos_l, Z.pow_pos_nonneg; [|apply Z.pow_nonneg]. | |
Qed. | |
Lemma pow_spec {x y}: pow x y = (x ^ y) mod modulus. | |
Proof. | |
case y; simpl; [| |reflexivity]. | |
- rewrite Zmod_1_l; [reflexivity|apply modulus_gt1]. | |
- intro p. unfold pow_pos, Z.pow_pos. | |
pattern p; apply Pos.peano_ind; [reflexivity|]. | |
clear p; intros p H; rewrite !Pos.iter_succ, H; unfold mul. | |
apply Zmult_mod_idemp_r. | |
Qed. | |
Lemma mul_succ x y: mul x (Z.succ y) = add x (mul x y). | |
Proof. | |
now unfold mul, add; rewrite Zplus_mod_idemp_r, Z.mul_succ_r, Z.add_comm. | |
Qed. | |
Lemma pow_succ {x y}: 0 <= y -> pow x (Z.succ y) = mul x (pow x y). | |
Proof. | |
intro H; rewrite !pow_spec, (Z.pow_succ_r _ _ H). | |
now unfold mul; rewrite Zmult_mod_idemp_r. | |
Qed. | |
Lemma pow_Zadd x y z: | |
0 <= y -> 0 <= z -> mul (pow x y) (pow x z) = pow x (y + z). | |
Proof. | |
unfold mul; rewrite !pow_spec; intros Hy Hz. | |
rewrite Zmult_mod_idemp_r, Zmult_mod_idemp_l, <-Z.pow_add_r; easy. | |
Qed. | |
Lemma Z2Nat_inj_pow x y: | |
0 <= x -> 0 <=y -> Z.to_nat (x ^ y) = (Z.to_nat x ^ Z.to_nat y)%nat. | |
Proof. | |
intros Hx Hy; pattern y; apply natlike_ind; auto. | |
clear y Hy; intros y Hy; simpl. | |
rewrite(Z2Nat.inj_succ y Hy),(Z.pow_succ_r _ y Hy); simpl. | |
rewrite(Z2Nat.inj_mul x _ Hx); simpl. intro E; rewrite E. | |
- reflexivity. | |
- apply Z.pow_nonneg; assumption. | |
Qed. | |
Lemma shiftl_spec {x y}: | |
0 <= x -> 0 <= y -> shiftl x y = Z.shiftl x y mod modulus. | |
Proof. | |
intros Hx Hy; unfold shiftl; rewrite Nat.shiftl_mul_pow2. | |
rewrite(Z.shiftl_mul_pow2 _ _ Hy). | |
assert(H2y: 0 <= 2 ^ y) by (apply Z.pow_nonneg; nia). | |
assert(Hxm: 0 <= x mod modulus) | |
by (apply Z_mod_lt; generalize modulus_gt1; nia). | |
replace 2%nat with (Z.to_nat 2) by easy. | |
rewrite <-Z2Nat_inj_pow by easy. | |
rewrite <-Z2Nat.inj_mul, Z2Nat.id, Zmult_mod_idemp_l; auto. | |
apply Z.mul_nonneg_nonneg; assumption. | |
Qed. | |
(* Axiom shiftl_mul_pow2 : forall x y, mul (pow 2 x) y = shiftl y x. *) | |
Lemma shiftl_mul_pow2 : forall x y, | |
0 <= x -> 0 <= y -> mul (pow 2 x) y = shiftl y x. | |
Proof. | |
intros x y Hx Hy. | |
unfold mul; rewrite pow_spec, Zmult_mod_idemp_l, (shiftl_spec Hy Hx). | |
now rewrite (Z.shiftl_mul_pow2 _ _ Hx), Z.mul_comm. | |
Qed. | |
(* Axiom pow_pow : forall x y z, pow (pow x y) z = pow x (mul y z). *) | |
Lemma pow_pow : forall x y z, | |
0 <= y -> 0 <= z -> pow (pow x y) z = pow x (y * z). | |
Proof. | |
intros x y z Hy Hz. | |
pattern z; apply natlike_ind; [now rewrite Z.mul_0_r| |assumption]. | |
clear z Hz; intros z Hz IH; rewrite(pow_succ Hz). | |
rewrite Z.mul_succ_r, IH, Z.add_comm, <-pow_Zadd; nia. | |
Qed. | |
Fixpoint Zseq x n := | |
match n with | |
| O => [] | |
| S n' => x :: Zseq (1 + x) n' | |
end. | |
Definition Zrange x y := | |
Zseq x (Z.to_nat (1 + y - x)). | |
Notation "[ x ;..; y ]" := (Zrange x y) (at level 0). | |
Lemma map2_ext : forall A B C (f g : A -> B -> C) xs ys, | |
(forall a b, f a b = g a b) -> | |
map2 f xs ys = map2 g xs ys. | |
Proof. | |
induction xs; simpl; intros. | |
auto. | |
destruct ys. | |
auto. | |
rewrite H. | |
f_equal. | |
auto. | |
Qed. | |
Instance map2_morphism A B C: | |
Proper (pointwise_relation _ (pointwise_relation _ eq) ==> @eq (list A) ==> @eq (list B) ==> @eq (list C)) (@map2 A B C). | |
Proof. | |
unfold Proper, respectful, pointwise_relation. | |
intros. | |
subst. | |
apply map2_ext; auto. | |
Qed. | |
Lemma map2_map : forall A A' B C (f : A -> B -> C) (g : A' -> A) xs ys, | |
map2 f (map g xs) ys = map2 (fun x y => f (g x) y) xs ys. | |
Proof. | |
induction xs; simpl; intros; | |
destruct ys; congruence. | |
Qed. | |
Definition fft_aux xs n := | |
let omegas := | |
map (fun k => pow (prt (n + 1)) k) [0 ;..; 2 ^ n] | |
in | |
map2 mul omegas xs. | |
Ltac setoid_autorewrite' rew := | |
lazymatch rew with | |
| (?rew',?lem',?lem) => setoid_rewrite lem + setoid_autorewrite' (rew',lem') | |
| (?lem',?lem) => setoid_rewrite lem || setoid_rewrite lem' | |
end. | |
Ltac setoid_autorewrite rew := repeat setoid_autorewrite' rew. | |
Ltac optimize' db rew := | |
eexists; intros; | |
match goal with | |
|- ?X = ?Y => | |
let P := fresh in | |
set (P := fun y => y = Y); | |
enough (P X) by auto | |
end; | |
repeat autounfold with db; | |
setoid_autorewrite rew; | |
reflexivity. | |
Hint Unfold fft_aux prt : my_db. | |
Definition fft_aux_sig : | |
{ f | forall xs n, fft_aux xs n = f xs n }. | |
Proof. | |
optimize' my_db (shiftl_mul_pow2,pow_pow,map2_map). | |
Defined. | |
Eval simpl in proj1_sig fft_aux_sig. | |
(* result : | |
[[[ | |
= fun (xs : list Z) (n : Z) => | |
map2 (fun x y : Z => shiftl y (shiftl x (N + 1 - (n + 1)))) [0;..; 2 ^ n] xs | |
: list Z -> Z -> list Z | |
]]] | |
*) | |
Function fft (xs : list Z) (n:nat) := | |
match n with | |
| O => xs | |
| S n' => | |
let (xs0,xs1) := divide_list xs in | |
let xs0' := fft xs0 n' in | |
let xs1' := proj1_sig fft_aux_sig (fft xs1 n') (Z.of_nat n') in | |
let xs' := map2 add xs0' xs1' in | |
let xs'' := map2 sub xs0' xs1' in | |
xs' ++ xs'' | |
end. | |
End Fft. |
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment