diff --git a/coq/Commands.v b/coq/Commands.v index 5999ba0d7db1ad75c2823dc3a39c83a3a37b9f3d..977733dff311f04eaa62205ac5ec054e26e573d8 100644 --- a/coq/Commands.v +++ b/coq/Commands.v @@ -7,8 +7,8 @@ Require Export Daisy.Infra.ExpressionAbbrevs Daisy.Infra.NatSet. (** Next define what a program is. - Currently no loops, only conditionals and assignments - Final return statement + Currently no loops, or conditionals. + Only assignments and return statement **) Inductive cmd (V:Type) :Type := Let: mType -> nat -> exp V -> cmd V -> cmd V @@ -43,7 +43,8 @@ Inductive sstep : cmd R -> env -> R -> cmd R -> env -> Prop := *) (** - Analogously define Big Step semantics for the Daisy language + Define big step semantics for the Daisy language, terminating on a "returned" + result value **) Inductive bstep : cmd R -> env -> R -> mType -> Prop := let_b m m1 x e s E v res: @@ -55,12 +56,19 @@ Inductive bstep : cmd R -> env -> R -> mType -> Prop := eval_exp E e v m -> bstep (Ret e) E v m. +(** + The free variables of a command are all used variables of expressions + without the let bound variables +**) Fixpoint freeVars V (f:cmd V) :NatSet.t := match f with - | Let _ x e1 g => NatSet.remove x (NatSet.union (Expressions.freeVars e1) (freeVars g)) - | Ret e => Expressions.freeVars e + | Let _ x e1 g => NatSet.remove x (NatSet.union (Expressions.usedVars e1) (freeVars g)) + | Ret e => Expressions.usedVars e end. +(** + The defined variables of a command are all let bound variables +**) Fixpoint definedVars V (f:cmd V) :NatSet.t := match f with | Let _ x _ g => NatSet.add x (definedVars g) diff --git a/coq/ErrorValidation.v b/coq/ErrorValidation.v index 8f963dea83493c6a0ea4370b09062aae1b32a499..18d28c4a1cffe113ce7a35ef023d18a3c2a16491 100644 --- a/coq/ErrorValidation.v +++ b/coq/ErrorValidation.v @@ -141,7 +141,7 @@ Proof. apply Rmult_le_compat_r. { apply mEps_geq_zero. } { rewrite <- maxAbs_impl_RmaxAbs. - apply contained_leq_maxAbs_val. + apply contained_leq_maxAbs. unfold contained; simpl. pose proof (validIntervalbounds_sound (Var Q x) A P (E:=fun y : nat => if y =? x then Some nR else E1 y) (vR:=nR) bounds_valid (fVars := (NatSet.add x fVars))) as valid_bounds_prf. rewrite absenv_var in valid_bounds_prf; simpl in valid_bounds_prf. diff --git a/coq/Expressions.v b/coq/Expressions.v index 28c445b0f4022e7bb551529b56a4256ab1a34234..e8422ff67f9099e6c58c539951d7fe8401552ff0 100644 --- a/coq/Expressions.v +++ b/coq/Expressions.v @@ -1,6 +1,5 @@ (** -Formalization of the base expression language for the daisy framework -Required in all files, since we will always reason about expressions. + Formalization of the base expression language for the daisy framework **) Require Import Coq.Reals.Reals Coq.micromega.Psatz Coq.QArith.QArith Coq.QArith.Qreals. Require Import Daisy.Infra.RealRationalProps. @@ -136,12 +135,9 @@ Definition perturb (r:R) (e:R) := (** Define expression evaluation relation parametric by an "error" epsilon. -This value will be used later to express float computations using a perturbation -of the real valued computation by (1 + delta), where |delta| <= machine epsilon. - -It is important that variables are not perturbed when loading from an environment. -This is the case, since loading a float value should not increase an additional error. -Unary negation is special! We do not have a new error here since IEE 754 gives us a sign bit +The result value expresses float computations according to the IEEE standard, +using a perturbation of the real valued computation by (1 + delta), where +|delta| <= machine epsilon. **) Inductive eval_exp (E:env) :(exp R) -> R -> mType -> Prop := | Var_load m m1 x v: @@ -172,17 +168,21 @@ Inductive eval_exp (E:env) :(exp R) -> R -> mType -> Prop := eval_exp E f1 v1 m1 -> eval_exp E (Downcast m f1) (perturb v1 delta) m. -Fixpoint freeVars (V:Type) (e:exp V) :NatSet.t := +(** + Define the set of "used" variables of an expression to be the set of variables + occuring in it +**) +Fixpoint usedVars (V:Type) (e:exp V) :NatSet.t := match e with | Var _ _ x => NatSet.singleton x - | Unop u e1 => freeVars e1 - | Binop b e1 e2 => NatSet.union (freeVars e1) (freeVars e2) - | Downcast _ e1 => freeVars e1 + | Unop u e1 => usedVars e1 + | Binop b e1 e2 => NatSet.union (usedVars e1) (usedVars e2) + | Downcast _ e1 => usedVars e1 | _ => NatSet.empty end. (** -If |delta| <= 0 then perturb v delta is exactly v. + If |delta| <= 0 then perturb v delta is exactly v. **) Lemma delta_0_deterministic (v:R) (delta:R): (Rabs delta <= 0)%R -> @@ -249,8 +249,7 @@ Qed. Helping lemma. Needed in soundness proof. For each evaluation of using an arbitrary epsilon, we can replace it by evaluating the subexpressions and then binding the result values to different -variables in the Eironment. -This relies on the property that variables are not perturbed as opposed to parameters +variables in the Environment. **) Lemma binary_unfolding b f1 f2 m E vF: eval_exp E (Binop b f1 f2) vF m -> @@ -302,7 +301,6 @@ Proof. (* auto. *) Qed. - (** Using the parametric expressions, define boolean expressions for conditionals **) @@ -311,7 +309,7 @@ Inductive bexp (V:Type) : Type := | less: exp V -> exp V -> bexp V. (** - Define evaluation of booleans for reals + Define evaluation of boolean expressions **) (* Inductive bval (E:env): (bexp R) -> Prop -> Prop := *) (* leq_eval (f1:exp R) (f2:exp R) (v1:R) (v2:R): *) diff --git a/coq/Infra/Abbrevs.v b/coq/Infra/Abbrevs.v index 89cc251c578b6766a6f1630b7002d092ae88b8b1..08c1f05436d16830c567a5ee47e97f5e6898add3 100644 --- a/coq/Infra/Abbrevs.v +++ b/coq/Infra/Abbrevs.v @@ -17,7 +17,6 @@ define them to automatically unfold upon simplification. **) Definition interval:Type := R * R. Definition err:Type := R. -Definition ann:Type := interval * err. Definition mkInterval (ivlo:R) (ivhi:R) := (ivlo,ivhi). Definition IVlo (intv:interval) := fst intv. Definition IVhi (intv:interval) := snd intv. @@ -37,9 +36,6 @@ Arguments mkIntv _ _/. Arguments ivlo _ /. Arguments ivhi _ /. -Ltac iv_assert iv name:= - assert (exists ivlo ivhi, iv = (ivlo, ivhi)) as name by (destruct iv; repeat eexists; auto). - (** Later we will argue about program preconditions. Define a precondition to be a function mapping numbers (resp. variables) to intervals. @@ -50,6 +46,10 @@ Definition precond :Type := nat -> intv. Abbreviation for the type of a variable environment, which should be a partial function **) Definition env := nat -> option (R * mType). + +(** + The empty environment must return NONE for every variable +**) Definition emptyEnv:env := fun _ => None. (** diff --git a/coq/Infra/Ltacs.v b/coq/Infra/Ltacs.v index c772c8388630b6526c5591c4b5c10e6b76b529ea..99dd3a84174ffa6bd39a4f2a476e74ca19ab7744 100644 --- a/coq/Infra/Ltacs.v +++ b/coq/Infra/Ltacs.v @@ -2,6 +2,10 @@ Require Import Coq.Bool.Bool Coq.Reals.Reals. Require Import Daisy.Infra.RealSimps Daisy.Infra.NatSet. + +Ltac iv_assert iv name:= + assert (exists ivlo ivhi, iv = (ivlo, ivhi)) as name by (destruct iv; repeat eexists; auto). + (** Automatic translation and destruction of conjuctinos with andb into Props **) Ltac andb_to_prop H := apply Is_true_eq_left in H; diff --git a/coq/ssaPrgs.v b/coq/ssaPrgs.v index 21c577cd9ad2bc438c4fb642a40bde47cc192319..7daf78c6d191373938e6ea1791bec792cb487f66 100644 --- a/coq/ssaPrgs.v +++ b/coq/ssaPrgs.v @@ -1,166 +1,111 @@ +(** + We define a pseudo SSA predicate. + The formalization is similar to the renamedApart property in the LVC framework by Schneider, Smolka and Hack + http://dblp.org/rec/conf/itp/SchneiderSH15 + Our predicate is not as fully fledged as theirs, but we especially borrow the idea of annotating + the program with the predicate with the set of free and defined variables +**) Require Import Coq.QArith.QArith Coq.QArith.Qreals Coq.Reals.Reals. Require Import Coq.micromega.Psatz. Require Import Daisy.Infra.RealRationalProps Daisy.Infra.Ltacs. Require Export Daisy.Commands. -Fixpoint validVars (V:Type) (f:exp V) Vs :bool := - match f with - | Const n => true - | Var _ _ v => NatSet.mem v Vs - | Unop o f1 => validVars f1 Vs - | Binop o f1 f2 => validVars f1 Vs && validVars f2 Vs - | Downcast _ f1 => validVars f1 Vs - end. - -Lemma validVars_subset_freeVars T (e:exp T) V: - validVars e V = true-> - NatSet.Subset (Expressions.freeVars e) V. -Proof. - revert V; induction e; simpl; intros V valid_V; try auto. - - rewrite NatSet.mem_spec in valid_V. - hnf. intros; rewrite NatSet.singleton_spec in *; subst; auto. - - hnf; intros a in_empty; inversion in_empty. - - andb_to_prop valid_V. - hnf; intros a in_union. - rewrite NatSet.union_spec in in_union. - destruct in_union as [in_e1 | in_e2]. - + specialize (IHe1 V L a in_e1); auto. - + specialize (IHe2 V R a in_e2); auto. -Qed. - -Lemma validVars_add V (f:exp V) Vs n: - validVars f Vs = true -> - validVars f (NatSet.add n Vs) = true. +Lemma validVars_add V (e:exp V) Vs n: + NatSet.Subset (usedVars e) Vs -> + NatSet.Subset (usedVars e) (NatSet.add n Vs). Proof. - induction f; try auto. - - intros valid_var. unfold validVars in *; simpl in *. - rewrite NatSet.mem_spec in *. + induction e; try auto. + - intros valid_subset. + hnf. intros a in_singleton. + specialize (valid_subset a in_singleton). rewrite NatSet.add_spec; right; auto. - intros vars_binop. simpl in *. - apply Is_true_eq_left in vars_binop. - apply Is_true_eq_true. - apply andb_prop_intro. - apply andb_prop_elim in vars_binop. - destruct vars_binop; split; apply Is_true_eq_left. - + apply IHf1. apply Is_true_eq_true; auto. - + apply IHf2. apply Is_true_eq_true; auto. -Qed. - -Inductive ssaPrg (V:Type): (cmd V) -> (NatSet.t) -> (NatSet.t) -> Prop := - ssaLet m x e s inVars Vterm: - validVars e inVars = true -> - NatSet.mem x inVars = false -> - ssaPrg s (NatSet.add x inVars) Vterm -> - ssaPrg (Let m x e s) inVars Vterm -|ssaRet e inVars Vterm: - validVars e inVars = true -> - NatSet.equal inVars (NatSet.remove 0%nat Vterm) = true-> - ssaPrg (Ret e) inVars Vterm. - - -Lemma validVars_valid_subset (V:Type) (e:exp V) inVars: - validVars e inVars = true -> - NatSet.Subset (Expressions.freeVars e) inVars. -Proof. - induction e; intros vars_valid; unfold validVars in *; simpl in *; try auto. - - rewrite NatSet.mem_spec in vars_valid. - hnf. intros; rewrite NatSet.singleton_spec in *; subst; auto. - - hnf; intros a in_empty; inversion in_empty. - - hnf; intros a in_union; rewrite NatSet.union_spec in in_union. - rewrite andb_true_iff in vars_valid. - destruct vars_valid. - destruct in_union as [in_e1 | in_e2]. - + apply IHe1; auto. - + apply IHe2; auto. + intros a in_empty. + inversion in_empty. + - simpl; intros in_vars. + intros a in_union. + rewrite NatSet.union_spec in in_union. + destruct in_union. + + apply IHe1; try auto. + hnf; intros. + apply in_vars. + rewrite NatSet.union_spec; auto. + + apply IHe2; try auto. + hnf; intros. + apply in_vars. + rewrite NatSet.union_spec; auto. Qed. -Lemma validVars_non_stuck (e:exp Q) inVars (E E':env): - validVars e inVars = true -> - E' = toREvalEnv E -> - (forall v, NatSet.In v (Expressions.freeVars e) -> exists r m, E' v = Some (r, m))%R -> - exists vR, eval_exp E' (toREval (toRExp e)) vR M0. +Lemma validVars_non_stuck (e:exp Q) inVars E: + NatSet.Subset (usedVars e) inVars -> + (forall v, NatSet.In v (usedVars e) -> + exists r, (toREvalEnv E) v = Some (r, M0))%R -> + exists vR, eval_exp (toREvalEnv E) (toREval (toRExp e)) vR M0. Proof. - revert inVars E; induction e; - intros inVars E vars_valid EnvR fVars_live. + revert inVars E; induction e; + intros inVars E vars_valid fVars_live. - simpl in *. assert (NatSet.In n (NatSet.singleton n)) as in_n by (rewrite NatSet.singleton_spec; auto). specialize (fVars_live n in_n). - destruct fVars_live as [vR [m1 E_def]]. - rewrite EnvR in E_def. - remember (E n) as en. - destruct en as [[r mr] | p]. - + unfold toREvalEnv in E_def. - rewrite <- Heqen in E_def. - inversion E_def. - subst. - exists vR. - pose proof (isMorePrecise_refl M0). - eapply (Var_load (toREvalEnv E) M0 n); eauto. - unfold toREvalEnv. - rewrite <- Heqen. - auto. - + unfold toREvalEnv in E_def. - rewrite <- Heqen in E_def. - inversion E_def. + destruct fVars_live as [vR E_def]. + exists vR; constructor; auto. - exists (perturb (Q2R v) 0); constructor. - simpl; rewrite Rabs_R0; rewrite Q2R0_is_0; lra. - - assert (exists vR, eval_exp E' (toREval (toRExp e)) vR M0) + simpl (meps M0); rewrite Rabs_R0; rewrite Q2R0_is_0; lra. + - assert (exists vR, eval_exp (toREvalEnv E) (toREval (toRExp e)) vR M0) as eval_e_def by (eapply IHe; eauto). destruct eval_e_def as [ve eval_e_def]. case_eq u; intros; subst. - + exists (evalUnop Neg ve); econstructor; eauto. - + exists (perturb (evalUnop Inv ve) 0); econstructor; eauto. - simpl. rewrite Q2R0_is_0. rewrite Rabs_R0; lra. - - andb_to_prop vars_valid; simpl in *. - assert (exists vR1, eval_exp E' (toREval (toRExp e1)) vR1 M0) as eval_e1_def. + + exists (evalUnop Neg ve); constructor; auto. + + exists (perturb (evalUnop Inv ve) 0); constructor; auto. + simpl (meps M0); rewrite Q2R0_is_0; rewrite Rabs_R0; lra. + - repeat rewrite NatSet.subset_spec in *; simpl in *. + assert (exists vR1, eval_exp (toREvalEnv E) (toREval (toRExp e1)) vR1 M0) as eval_e1_def. + eapply IHe1; eauto. - intros. - destruct (fVars_live v) as [vR E_def]; try eauto. - apply NatSet.union_spec; auto. - + assert (exists vR2, eval_exp E' (toREval (toRExp e2)) vR2 M0) as eval_e2_def. - * eapply IHe2; eauto. - intros. + * hnf; intros. + apply vars_valid. + rewrite NatSet.union_spec; auto. + * intros. destruct (fVars_live v) as [vR E_def]; try eauto. apply NatSet.union_spec; auto. + + assert (exists vR2, eval_exp (toREvalEnv E) (toREval (toRExp e2)) vR2 M0) as eval_e2_def. + * eapply IHe2; eauto. + { hnf; intros. + apply vars_valid. + rewrite NatSet.union_spec; auto. } + { intros. + destruct (fVars_live v) as [vR E_def]; try eauto. + apply NatSet.union_spec; auto. } * destruct eval_e1_def as [vR1 eval_e1_def]; destruct eval_e2_def as [vR2 eval_e2_def]. exists (perturb (evalBinop b vR1 vR2) 0); econstructor; eauto. auto. simpl. rewrite Q2R0_is_0. rewrite Rabs_R0; lra. - - assert (exists vR, eval_exp E' (toREval (toRExp e)) vR M0) as eval_r_def by (eapply IHe; eauto). + - assert (exists vR, eval_exp (toREvalEnv E) (toREval (toRExp e)) vR M0) as eval_r_def by (eapply IHe; eauto). destruct eval_r_def as [vr eval_r_def]. exists vr. simpl toREval. auto. Qed. -Lemma validVars_equal_set V (e:exp V) vars vars': - NatSet.Equal vars vars' -> - validVars e vars = true -> - validVars e vars' = true. -Proof. - revert vars vars'. - induction e; intros vars vars' eq_set valid_vars; simpl in *; auto. - - rewrite NatSet.mem_spec in *. - rewrite <- eq_set; auto. - - eapply IHe; eauto. - - apply Is_true_eq_true. - apply andb_prop_intro. - andb_to_prop valid_vars. - split; apply Is_true_eq_left. - + eapply IHe1; eauto. - + eapply IHe2; eauto. - - eapply IHe; eauto. -Qed. +Inductive ssaPrg (V:Type): (cmd V) -> (NatSet.t) -> (NatSet.t) -> Prop := + ssaLet m x e s inVars Vterm: + NatSet.Subset (usedVars e) inVars-> + NatSet.mem x inVars = false -> + ssaPrg s (NatSet.add x inVars) Vterm -> + ssaPrg (Let m x e s) inVars Vterm +|ssaRet e inVars Vterm: + NatSet.Subset (usedVars e) inVars -> + NatSet.Equal inVars Vterm -> + ssaPrg (Ret e) inVars Vterm. + Lemma ssa_subset_freeVars V (f:cmd V) inVars outVars: ssaPrg f inVars outVars -> NatSet.Subset (Commands.freeVars f) inVars. Proof. intros ssa_f; induction ssa_f. - - simpl in *. apply validVars_subset_freeVars in H. - hnf; intros a in_fVars. + - simpl in *. hnf; intros a in_fVars. rewrite NatSet.remove_spec, NatSet.union_spec in in_fVars. destruct in_fVars as [in_union not_eq]. destruct in_union; try auto. @@ -169,7 +114,6 @@ Proof. destruct IHssa_f; subst; try auto. exfalso; apply not_eq; auto. - hnf; intros. - apply validVars_subset_freeVars in H. simpl in H1. apply H; auto. Qed. @@ -184,8 +128,7 @@ Proof. revert set_eq; revert inVars'. induction ssa_f. - constructor. - + eapply validVars_equal_set; eauto. - symmetry; auto. + + rewrite set_eq; auto. + case_eq (NatSet.mem x inVars'); intros case_mem; try auto. rewrite NatSet.mem_spec in case_mem. rewrite set_eq in case_mem. @@ -194,21 +137,16 @@ Proof. + apply IHssa_f; auto. apply NatSetProps.Dec.F.add_m; auto. - constructor. - + eapply validVars_equal_set; eauto. - symmetry; auto. - + rewrite NatSet.equal_spec in *. - hnf. intros a; split. - * intros in_primed. - rewrite <- H0. rewrite <- set_eq. auto. - * intros in_rem. rewrite set_eq. rewrite H0; auto. + + rewrite set_eq; auto. + + rewrite set_eq; auto. Qed. Fixpoint validSSA (f:cmd Q) (inVars:NatSet.t) := match f with |Let m x e g => - andb (andb (negb (NatSet.mem x inVars)) (validVars e inVars)) (validSSA g (NatSet.add x inVars)) - |Ret e => validVars e inVars && (negb (NatSet.mem 0%nat inVars)) + andb (andb (negb (NatSet.mem x inVars)) (NatSet.subset (usedVars e) inVars)) (validSSA g (NatSet.add x inVars)) + |Ret e => NatSet.subset (usedVars e) inVars end. Lemma validSSA_sound f inVars: @@ -223,32 +161,14 @@ Proof. destruct IHf as [outVars IHf]. exists outVars. constructor; eauto. - rewrite negb_true_iff in L0. auto. + + rewrite <- NatSet.subset_spec; auto. + + rewrite negb_true_iff in L0. auto. - intros inVars validSSA_ret. simpl in *. - exists (NatSet.add 0%nat inVars). - andb_to_prop validSSA_ret. + exists inVars. constructor; auto. - rewrite negb_true_iff in R. - hnf in R. - rewrite NatSet.equal_spec. - hnf. - intros a. - rewrite NatSet.remove_spec, NatSet.add_spec. - split. - + intros in_inVars. - case_eq (a =? 0%nat). - * intros a_zero. - rewrite Nat.eqb_eq in a_zero. - rewrite a_zero in in_inVars. - rewrite <- NatSet.mem_spec in in_inVars. - rewrite in_inVars in R. - inversion R. - * intros a_neq_zero. apply beq_nat_false in a_neq_zero. - split; auto. - + intros in_add_rem. - destruct in_add_rem as [ [a_zero | a_inVars] a_neq_zero]; try auto. - exfalso; eauto. + + rewrite <- NatSet.subset_spec; auto. + + hnf; intros; split; auto. Qed. Lemma ssa_shadowing_free m x y v v' e f inVars outVars E: @@ -341,7 +261,7 @@ revert E1 E2 vR. Qed. Lemma dummy_bind_ok e v x v' inVars E: - validVars e inVars = true -> + NatSet.Subset (usedVars e) inVars -> NatSet.mem x inVars = false -> eval_exp E (toREval (toRExp e)) v M0 -> eval_exp (updEnv x M0 v' E) (toREval (toRExp e)) v M0. @@ -355,25 +275,32 @@ Proof. intros n_eq_x. rewrite Nat.eqb_eq in n_eq_x. subst; simpl in *. - rewrite x_not_free in valid_vars; inversion valid_vars. + hnf in valid_vars. + assert (NatSet.mem x (NatSet.singleton x) = true) + as in_singleton by (rewrite NatSet.mem_spec, NatSet.singleton_spec; auto). + rewrite NatSet.mem_spec in *. + specialize (valid_vars x in_singleton). + rewrite <- NatSet.mem_spec in valid_vars. + rewrite valid_vars in *; congruence. + econstructor. auto. rewrite H; auto. - inversion eval_e; subst; constructor; auto. - inversion eval_e; subst; econstructor; eauto. - simpl in valid_vars. - apply Is_true_eq_left in valid_vars. - apply andb_prop_elim in valid_vars. - destruct valid_vars. inversion eval_e; subst; econstructor; eauto; assert (M0 = M0) as M00 by auto; - pose proof (ifM0isJoin_l M0 m1 m2 M00 H4); - pose proof (ifM0isJoin_r M0 m1 m2 M00 H4); + pose proof (ifM0isJoin_l M0 m1 m2 M00 H2); + pose proof (ifM0isJoin_r M0 m1 m2 M00 H2); subst. + eapply IHe1; eauto. - apply Is_true_eq_true; auto. + hnf; intros a in_e1. + apply valid_vars; + rewrite NatSet.union_spec; auto. + eapply IHe2; eauto. - apply Is_true_eq_true; auto. + hnf; intros a in_e2. + apply valid_vars; + rewrite NatSet.union_spec; auto. - apply (IHe v1 x v2 inVars E); auto. Qed. @@ -507,7 +434,7 @@ Admitted. Lemma stepwise_substitution x e v f E vR inVars outVars: ssaPrg (toREvalCmd (toRCmd f)) inVars outVars -> NatSet.In x inVars -> - validVars e inVars = true -> + NatSet.Subset (usedVars e) inVars -> eval_exp (toREvalEnv E) (toREval (toRExp e)) v M0 -> bstep (toREvalCmd (toRCmd f)) (updEnv x M0 v (toREvalEnv E)) vR M0 <-> bstep (toREvalCmd (toRCmd (map_subst f x e))) (toREvalEnv E) vR M0. diff --git a/hol4/abbrevsScript.sml b/hol4/AbbrevsScript.sml similarity index 56% rename from hol4/abbrevsScript.sml rename to hol4/AbbrevsScript.sml index 00ea9353b4154d6f0a25d20af3e66eafae7bad22..ce477c1e48fa47161342e4ec937909bf6f638beb 100644 --- a/hol4/abbrevsScript.sml +++ b/hol4/AbbrevsScript.sml @@ -1,7 +1,10 @@ +(** + This file contains some type abbreviations, to ease writing. + **) open preamble -open realTheory realLib +open realTheory realLib sptreeTheory -val _ = new_theory "abbrevs"; +val _ = new_theory "Abbrevs"; (** For the moment we need only one interval type in HOL, since we do not have the problem of computability as we have it in Coq @@ -10,8 +13,6 @@ val _ = type_abbrev("interval", ``:real#real``); val IVlo_def = Define `IVlo (iv:interval) = FST iv`; val IVhi_def = Define `IVhi (iv:interval) = SND iv`; -val _ = type_abbrev("ann", ``:interval#real``); - (** Later we will argue about program preconditions. Define a precondition to be a function mapping numbers (resp. variables) to intervals. @@ -19,15 +20,25 @@ Define a precondition to be a function mapping numbers (resp. variables) to inte val _ = type_abbrev ("precond", ``:num->interval``); (** - Abbreviations for the environment types + Abbreviation for the type of a variable environment, which should be a partial function **) val _ = type_abbrev("env",``:num->real option``); (** - Define environment update function as abbreviation. + The empty environment must return NONE for every variable +**) +val emptyEnv_def = Define ` + emptyEnv (x:num) = NONE`; + +(** + Define environment update function as abbreviation, for variable environments **) val updEnv_def = Define ` updEnv (x:num) (v:real) (E:env) (y:num) :real option = if y = x then SOME v else E y`; +(** + Abbreviation for insertion into "num_set" type +**) +val Insert_def = Define `Insert x V = insert x () V`; val - = export_theory(); diff --git a/hol4/CertificateCheckerScript.sml b/hol4/CertificateCheckerScript.sml index 3d5beefc1d7b519feed86a2abba1868a29690ee0..ef11474e7ab660735c2d2e034d0e9e2255ded562 100644 --- a/hol4/CertificateCheckerScript.sml +++ b/hol4/CertificateCheckerScript.sml @@ -6,7 +6,7 @@ **) open preamble open simpLib realTheory realLib RealArith -open abbrevsTheory expressionsTheory RealSimpsTheory +open AbbrevsTheory ExpressionsTheory RealSimpsTheory open ExpressionAbbrevsTheory ErrorBoundsTheory IntervalArithTheory open IntervalValidationTheory ErrorValidationTheory diff --git a/hol4/CommandsScript.sml b/hol4/CommandsScript.sml index c477f929cd9e3e9fc046ccf1c80bb86237e6ed6c..e15f63c9a331374251341a1cb355543bd9821cc9 100644 --- a/hol4/CommandsScript.sml +++ b/hol4/CommandsScript.sml @@ -3,7 +3,7 @@ **) open preamble open simpLib realTheory realLib RealArith -open abbrevsTheory ExpressionAbbrevsTheory +open AbbrevsTheory ExpressionAbbrevsTheory val _ = new_theory "Commands"; diff --git a/hol4/EnvironmentsScript.sml b/hol4/EnvironmentsScript.sml index 5f2ae331b211454c861e166b0b06f4473ada177c..bafe0ca87772402e4e3ec6082663528b220c4b3d 100644 --- a/hol4/EnvironmentsScript.sml +++ b/hol4/EnvironmentsScript.sml @@ -1,6 +1,6 @@ open preamble open simpLib realTheory realLib RealArith -open abbrevsTheory ExpressionAbbrevsTheory CommandsTheory +open AbbrevsTheory ExpressionAbbrevsTheory CommandsTheory val _ = new_theory "Environments"; diff --git a/hol4/ErrorBoundsScript.sml b/hol4/ErrorBoundsScript.sml index 0dd174f62aa2537fbf9e8c46796a77f2dd17b5a5..62d4db05438afd1bf4ac1d64de5c732ee03fa11d 100644 --- a/hol4/ErrorBoundsScript.sml +++ b/hol4/ErrorBoundsScript.sml @@ -1,11 +1,11 @@ (** -Proofs of general bounds on the error of arithmetic expressions. +Proofs of general bounds on the error of arithmetic Expressions. This shortens soundness proofs later. Bounds are explained in section 5, Deriving Computable Error Bounds **) open preamble open simpLib realTheory realLib RealArith -open abbrevsTheory expressionsTheory RealSimpsTheory +open AbbrevsTheory ExpressionsTheory RealSimpsTheory open ExpressionAbbrevsTheory EnvironmentsTheory val _ = new_theory "ErrorBounds"; diff --git a/hol4/ErrorValidationScript.sml b/hol4/ErrorValidationScript.sml index baff42ededcb58657882c3a85a80b6d1d879c2db..8fbc9e7530c82a495da0d6b75da068db6677832f 100644 --- a/hol4/ErrorValidationScript.sml +++ b/hol4/ErrorValidationScript.sml @@ -8,7 +8,7 @@ **) open preamble open simpLib realTheory realLib RealArith -open abbrevsTheory expressionsTheory RealSimpsTheory +open AbbrevsTheory ExpressionsTheory RealSimpsTheory open ExpressionAbbrevsTheory ErrorBoundsTheory IntervalArithTheory open IntervalValidationTheory EnvironmentsTheory CommandsTheory ssaPrgsTheory diff --git a/hol4/ExpressionAbbrevsScript.sml b/hol4/ExpressionAbbrevsScript.sml index 3605b0f4a3b576156e8610d69a5fa811641d5292..cf85405dd998e2f5e9ddb3c2977f13f73e1e0769 100644 --- a/hol4/ExpressionAbbrevsScript.sml +++ b/hol4/ExpressionAbbrevsScript.sml @@ -1,6 +1,6 @@ open preamble -open expressionsTheory +open ExpressionsTheory val _ = new_theory "ExpressionAbbrevs" diff --git a/hol4/expressionsScript.sml b/hol4/ExpressionsScript.sml similarity index 54% rename from hol4/expressionsScript.sml rename to hol4/ExpressionsScript.sml index 3f4746df3780b90c2e146333e46a8cbfe8cbfe28..fbd5b287a3d18bad1739ba5f108e95470d88066b 100644 --- a/hol4/expressionsScript.sml +++ b/hol4/ExpressionsScript.sml @@ -1,9 +1,13 @@ -open preamble -open realTheory realLib -open abbrevsTheory +(** + Formalization of the base expression language for the daisy framework + **) +open preamble miscTheory +open realTheory realLib sptreeTheory +open AbbrevsTheory + val _ = ParseExtras.temp_tight_equality() -val _ = new_theory "expressions"; +val _ = new_theory "Expressions"; (** Expressions will use binary operators. @@ -36,14 +40,10 @@ Errors are added in the expression evaluation level later val evalUnop_def = Define ` evalUnop Neg v = - v /\ evalUnop Inv v = inv(v:real)` + (** - Define expressions parametric over some value type 'v. Will ease - reasoning about different instantiations later. Note that we - differentiate between whether we use a variable from the environment - or a parameter. Parameters do not have error bounds in the - invariants, so they must be perturbed, but variables from the - program will be perturbed upon binding, so we do not need to perturb - them. + Define Expressions parametric over some value type 'v. + Will ease reasoning about different instantiations later. **) val _ = Datatype ` exp = Var num @@ -58,30 +58,34 @@ val perturb_def = Define ` perturb (r:real) (e:real) = r * (1 + e)` (** - Define expression evaluation relation parametric by an "error" - delta. This value will be used later to express float computations - using a perturbation of the real valued computation by (1 + d) +Define expression evaluation relation parametric by an "error" epsilon. +The result value expresses float computations according to the IEEE standard, +using a perturbation of the real valued computation by (1 + delta), where +|delta| <= machine epsilon. **) val (eval_exp_rules, eval_exp_ind, eval_exp_cases) = Hol_reln ` - (!eps E x v. - E x = SOME v ==> - eval_exp eps E (Var x) v) /\ - (!eps E n delta. - abs delta <= eps ==> - eval_exp eps E (Const n) (perturb n delta)) /\ - (!eps E f1 v1. - eval_exp eps E f1 v1 ==> - eval_exp eps E (Unop Neg f1) (evalUnop Neg v1)) /\ - (!eps E f1 v1 delta. - abs delta <= eps /\ - eval_exp eps E f1 v1 ==> - eval_exp eps E (Unop Inv f1) (perturb (evalUnop Inv v1) delta)) /\ - (!eps E b f1 f2 v1 v2 delta. - abs delta <= eps /\ - eval_exp eps E f1 v1 /\ - eval_exp eps E f2 v2 ==> - eval_exp eps E (Binop b f1 f2) (perturb (evalBinop b v1 v2) delta))`; + (!eps E x v. + E x = SOME v ==> + eval_exp eps E (Var x) v) /\ + (!eps E n delta. + abs delta <= eps ==> + eval_exp eps E (Const n) (perturb n delta)) /\ + (!eps E f1 v1. + eval_exp eps E f1 v1 ==> + eval_exp eps E (Unop Neg f1) (evalUnop Neg v1)) /\ + (!eps E f1 v1 delta. + abs delta <= eps /\ + eval_exp eps E f1 v1 ==> + eval_exp eps E (Unop Inv f1) (perturb (evalUnop Inv v1) delta)) /\ + (!eps E b f1 f2 v1 v2 delta. + abs delta <= eps /\ + eval_exp eps E f1 v1 /\ + eval_exp eps E f2 v2 ==> + eval_exp eps E (Binop b f1 f2) (perturb (evalBinop b v1 v2) delta))`; +(** + Generate a better case lemma +**) val eval_exp_cases = map (GEN_ALL o SIMP_CONV (srw_ss()) [Once eval_exp_cases]) [``eval_exp eps E (Var v) res``, @@ -90,10 +94,28 @@ val eval_exp_cases = ``eval_exp eps E (Binop n e1 e2) res``] |> LIST_CONJ |> curry save_thm "eval_exp_cases"; +(** + Define the set of "used" variables of an expression to be the set of variables + occuring in it +**) +val usedVars_def = Define ` + usedVars (e: 'a exp) :num_set = + case e of + | Var x => insert x () (LN) + | Unop u e1 => usedVars e1 + | Binop b e1 e2 => union (usedVars e1) (usedVars e2) + | _ => LN`; + +(** + If |delta| <= 0 then perturb v delta is exactly v. +**) val delta_0_deterministic = store_thm("delta_0_deterministic", ``!(v:real) (delta:real). abs delta <= 0 ==> perturb v delta = v``, fs [perturb_def,ABS_BOUNDS,REAL_LE_ANTISYM]); +(** +Evaluation with 0 as machine epsilon is deterministic +**) val meps_0_deterministic = store_thm("meps_0_deterministic", ``!(e: real exp) v1:real v2:real E. eval_exp (&0) E e v1 /\ eval_exp (&0) E e v2 ==> v1 = v2``, @@ -104,22 +126,33 @@ val meps_0_deterministic = store_thm("meps_0_deterministic", (** Helping lemma. Needed in soundness proof. For each evaluation of using an arbitrary epsilon, we can replace it by -evaluating the subexpressions and then binding the result values to different +evaluating the subExpressions and then binding the result values to different variables in the Eironment. -This relies on the property that variables are not perturbed as opposed to parameters **) val binary_unfolding = store_thm("binary_unfolding", -``!(b:binop) (e1:(real)exp) (e2:(real)exp) (eps:real) E (v:real). + ``!(b:binop) (e1:(real)exp) (e2:(real)exp) (eps:real) E (v:real). (eval_exp eps E (Binop b e1 e2) v <=> - (?(v1:real) (v2:real). - eval_exp eps E e1 v1 /\ - eval_exp eps E e2 v2 /\ - eval_exp eps (updEnv 2 v2 (updEnv 1 v1 E)) (Binop b (Var 1) (Var 2)) v))``, + (?(v1:real) (v2:real). + eval_exp eps E e1 v1 /\ + eval_exp eps E e2 v2 /\ + eval_exp eps (updEnv 2 v2 (updEnv 1 v1 emptyEnv)) (Binop b (Var 1) (Var 2)) v))``, + fs [updEnv_def, eval_exp_cases,APPLY_UPDATE_THM,PULL_EXISTS] + \\ metis_tac []); + +(** + Analogous lemma for unary expressions +**) +val unary_unfolding = store_thm("unary_unfolding", + ``!(u:unop) (e1:(real)exp) (eps:real) E (v:real). + (eval_exp eps E (Unop Inv e1) v <=> + (?(v1:real). + eval_exp eps E e1 v1 /\ + eval_exp eps (updEnv 1 v1 emptyEnv) (Unop Inv (Var 1)) v))``, fs [updEnv_def, eval_exp_cases,APPLY_UPDATE_THM,PULL_EXISTS] \\ metis_tac []); (* - Using the parametric expressions, define boolean expressions for conditionals + Using the parametric Expressions, define boolean Expressions for conditionals *) val _ = Datatype ` bexp = Leq ('v exp) ('v exp) @@ -144,7 +177,9 @@ val bval_exp_cases = ``bval eps E (Less e1 e2) res``] |> LIST_CONJ |> curry save_thm "bval_exp_cases"; -(* Simplify arithmetic later by making > >= only abbreviations *) +(** + Simplify arithmetic later by making > >= only abbreviations +**) val _ = overload_on("Gr",``\(e1:('v)exp). \(e2:('v)exp). Less e2 e1``); val _ = overload_on("Greq",``\(e1:('v)exp). \(e2:('v)exp). Leq e2 e1``); diff --git a/hol4/Holmakefile b/hol4/Holmakefile index a0d3397708b3ec1cf3ef6eb5928da29b6110d7e7..95bbc4f814b07d749df5a27bfbe8e3aad12b3992 100644 --- a/hol4/Holmakefile +++ b/hol4/Holmakefile @@ -1,4 +1,4 @@ -INCLUDES = $(HOLDIR)/examples/machine-code/hoare-triple common +INCLUDES = $(HOLDIR)/examples/machine-code/hoare-triple Infra OPTIONS = QUIT_ON_FAILURE ifdef POLY @@ -18,7 +18,7 @@ BARE_THYS = BasicProvers Defn HolKernel Parse Tactic monadsyntax \ optionTheory pairLib pairTheory pred_setTheory \ quantHeuristicsLib relationTheory res_quanTheory rich_listTheory \ sortingTheory sptreeTheory stringTheory sumTheory wordsTheory \ - simpLib realTheory realLib RealArith common/preamble + simpLib realTheory realLib RealArith Infra/preamble DEPS = $(patsubst %,%.uo,$(BARE_THYS1)) $(PARENTHEAP) diff --git a/hol4/common/Holmakefile b/hol4/Infra/Holmakefile similarity index 100% rename from hol4/common/Holmakefile rename to hol4/Infra/Holmakefile diff --git a/hol4/common/miscScript.sml b/hol4/Infra/miscScript.sml similarity index 100% rename from hol4/common/miscScript.sml rename to hol4/Infra/miscScript.sml diff --git a/hol4/common/preamble.sml b/hol4/Infra/preamble.sml similarity index 100% rename from hol4/common/preamble.sml rename to hol4/Infra/preamble.sml diff --git a/hol4/IntervalArithScript.sml b/hol4/IntervalArithScript.sml index 02b101675424dcb614278f4323a963205cff3f6c..239f98e73e51c1431c29acb9572e13497de86fe3 100644 --- a/hol4/IntervalArithScript.sml +++ b/hol4/IntervalArithScript.sml @@ -4,7 +4,7 @@ **) open preamble open realTheory realLib RealArith -open abbrevsTheory expressionsTheory RealSimpsTheory; +open AbbrevsTheory ExpressionsTheory RealSimpsTheory; val _ = new_theory "IntervalArith"; diff --git a/hol4/IntervalValidationScript.sml b/hol4/IntervalValidationScript.sml index fc311b71e097ad1eda1ed8a3263d27f15829e5cc..edce502eb446fec9c31172b4871308d457268051 100644 --- a/hol4/IntervalValidationScript.sml +++ b/hol4/IntervalValidationScript.sml @@ -7,7 +7,7 @@ **) open preamble open simpLib realTheory realLib RealArith -open abbrevsTheory expressionsTheory RealSimpsTheory +open AbbrevsTheory ExpressionsTheory RealSimpsTheory open ExpressionAbbrevsTheory IntervalArithTheory CommandsTheory ssaPrgsTheory val _ = new_theory "IntervalValidation"; diff --git a/src/main/scala/daisy/Context.scala b/src/main/scala/daisy/Context.scala index 73b0fd8f71607dc275cbbd2fef0aa802e71f1b10..cb0edc41a9b194f0283a5629d2e86d36337f6972 100644 --- a/src/main/scala/daisy/Context.scala +++ b/src/main/scala/daisy/Context.scala @@ -12,6 +12,8 @@ import scala.collection.immutable.Seq import scala.reflect.ClassTag import lang.Identifiers._ +import lang.Trees._ + import utils.{Interval, PartialInterval, Rational} case class Context( @@ -28,15 +30,22 @@ case class Context( // but don't want to pollute the nice and clean trees. // If these get too many, move to their own "Summary". // indexed by FunDef.id - inputRanges: Map[Identifier, Map[Identifier, Interval]] = Map(), - inputErrors: Map[Identifier, Map[Identifier, Rational]] = Map(), + specInputRanges: Map[Identifier, Map[Identifier, Interval]] = Map(), + specInputErrors: Map[Identifier, Map[Identifier, Rational]] = Map(), // for now we only support a single result value, i.e. no tuples // this map is indexed by fnc.id -> potentially partial interval bound of result // and similar for the errors - resultRangeBounds: Map[Identifier, PartialInterval] = Map(), - resultErrorBounds: Map[Identifier, Rational] = Map() - //requiredOutputRanges: Map[Identifier, Map[Identifier, PartialInterval]] = Map(), - //requiredOutputErrors: Map[Identifier, Map[Identifier, Rational]] = Map() + specResultRangeBounds: Map[Identifier, PartialInterval] = Map(), + specResultErrorBounds: Map[Identifier, Rational] = Map(), + + // the intermediate analysis results, since we rely on having them at hand when + // encoding the analysis result + intermediateAbsoluteErrors: Map[Identifier, Map[Expr, Rational]]= Map(), + intermediateRanges: Map[Identifier, Map[Expr, Interval]] = Map(), + + // the analysed/computed roundoff errors for each function + resultAbsoluteErrors: Map[Identifier, Rational] = Map(), + resultRealRanges: Map[Identifier, Interval] = Map() ) { // on the first creation of a context, we also update the context variable diff --git a/src/main/scala/daisy/InfoPhase.scala b/src/main/scala/daisy/InfoPhase.scala index bb4d2eef2383ebe1a0c275126053a872d4cf8137..37a977dc119279c2e02a314b80482f4f58702e91 100644 --- a/src/main/scala/daisy/InfoPhase.scala +++ b/src/main/scala/daisy/InfoPhase.scala @@ -50,23 +50,16 @@ object InfoPhase extends DaisyPhase { reporter.info(fnc.id) - val infoString: String = getLastExpression(fnc.body.get) match { - case x: NumAnnotation if x.hasError => - val absError = x.absError - val range = x.interval - val relError: Double = if (range.xlo <= zero && zero <= range.xhi) { - Double.NaN - } else { - max(abs(absError / range.xlo), abs(absError / range.xhi)).toDouble - } - - s"abs-error: ${x.absError}, range: ${x.interval},\nrel-error: ${relError}" - - case x: NumAnnotation => - "final error: - " - - case _ => "" - } + val absError = ctx.resultAbsoluteErrors(fnc.id) + val range = ctx.resultRealRanges(fnc.id) + val relError: Double = if (range.xlo <= zero && zero <= range.xhi) { + Double.NaN + } else { + max(abs(absError / range.xlo), abs(absError / range.xhi)).toDouble + } + + val infoString: String = + s"abs-error: ${absError}, real range: ${range},\nrel-error: ${relError}" reporter.info(infoString) if (out != null) { diff --git a/src/main/scala/daisy/analysis/DynamicPhase.scala b/src/main/scala/daisy/analysis/DynamicPhase.scala index 52cd0099a4fbb89c01ff1fa2eb312037b6334919..5161a7975528e881bb95a4e8d3226850ef3a2b85 100644 --- a/src/main/scala/daisy/analysis/DynamicPhase.scala +++ b/src/main/scala/daisy/analysis/DynamicPhase.scala @@ -106,7 +106,7 @@ object DynamicPhase extends DaisyPhase { val body = fnc.body.get reporter.info("evaluating " + id + "...") - val inputRanges: Map[Identifier, Interval] = ctx.inputRanges(id).map({ + val inputRanges: Map[Identifier, Interval] = ctx.specInputRanges(id).map({ case (id, i) => (id, Interval(i.mid - inputRangeFactor * i.radius, i.mid + inputRangeFactor * i.radius)) diff --git a/src/main/scala/daisy/analysis/RangeErrorPhase.scala b/src/main/scala/daisy/analysis/RangeErrorPhase.scala index b66c99d3f4adf05f44af951eff8cc6957313ee96..bdc565db33bb3ec01da4ca92c8bc3edb6406fc52 100644 --- a/src/main/scala/daisy/analysis/RangeErrorPhase.scala +++ b/src/main/scala/daisy/analysis/RangeErrorPhase.scala @@ -16,7 +16,7 @@ import lang.Identifiers._ Prerequisites: - SpecsProcessingPhase */ -object RangeErrorPhase extends DaisyPhase with ErrorFunctions { +object RangeErrorPhase extends DaisyPhase with RoundoffEvaluators with IntervalSubdivision { override val name = "range-error phase" override val description = "Computes ranges and absolute errors" @@ -37,10 +37,6 @@ object RangeErrorPhase extends DaisyPhase with ErrorFunctions { var reporter: Reporter = null - // Var for error functions - // trackRoundoffErrs and uniformPrecision assigned below in run - attachToTree = true - override def run(ctx: Context, prg: Program): (Context, Program) = { reporter = ctx.reporter reporter.info(s"\nStarting $name") @@ -56,9 +52,8 @@ object RangeErrorPhase extends DaisyPhase with ErrorFunctions { var errorMethod = "interval" var trackInitialErrs = true + var trackRoundoffErrs = true - // setting trait variables - trackRoundoffErrs = true var uniformPrecision: Precision = Float64 // process relevant options @@ -97,11 +92,16 @@ object RangeErrorPhase extends DaisyPhase with ErrorFunctions { val fncsToConsider: Seq[String] = functionsToConsider(ctx, prg) - for (fnc <- prg.defs) - if (!fnc.precondition.isEmpty && !fnc.body.isEmpty && fncsToConsider.contains(fnc.id.toString)){ + var roundoffErrorMap: Map[Identifier, Map[Expr, Rational]] = Map() + var rangeResMap: Map[Identifier, Map[Expr, Interval]] = Map() + + val res: Map[Identifier, (Rational, Interval)] = prg.defs.filter(fnc => + !fnc.precondition.isEmpty && + !fnc.body.isEmpty && + fncsToConsider.contains(fnc.id.toString)).map(fnc => { reporter.info("analyzing fnc: " + fnc.id) - val inputValMap: Map[Identifier, Interval] = ctx.inputRanges(fnc.id) + val inputValMap: Map[Identifier, Interval] = ctx.specInputRanges(fnc.id) // If we track both input and roundoff errors, then we pre-compute // the roundoff errors for those variables that do not have a user-defined @@ -109,14 +109,14 @@ object RangeErrorPhase extends DaisyPhase with ErrorFunctions { val inputErrorMap: Map[Identifier, Rational] = if (trackInitialErrs && trackRoundoffErrs){ - val inputErrs = ctx.inputErrors(fnc.id) + val inputErrs = ctx.specInputErrors(fnc.id) val allIDs = fnc.params.map(_.id).toSet val missingIDs = allIDs -- inputErrs.keySet inputErrs ++ missingIDs.map( id => (id -> uniformPrecision.absRoundoff(inputValMap(id)))) } else if(trackInitialErrs) { - val inputErrs = ctx.inputErrors(fnc.id) + val inputErrs = ctx.specInputErrors(fnc.id) val allIDs = fnc.params.map(_.id).toSet val missingIDs = allIDs -- inputErrs.keySet inputErrs ++ missingIDs.map( id => (id -> zero)) @@ -133,37 +133,46 @@ object RangeErrorPhase extends DaisyPhase with ErrorFunctions { } - (rangeMethod, errorMethod) match { + val (resError: Rational, resRange: Interval, rangeMap:Map[Expr, Interval], errorMap:Map[Expr, Interval]) = (rangeMethod, errorMethod) match { case ("interval", "interval") => - errorIntervalInterval(fnc.body.get, inputValMap, inputErrorMap) + uniformRoundoff_IA_IA(fnc.body.get, inputValMap, inputErrorMap, uniformPrecision, trackRoundoffErrs) case ("interval", "affine") => - errorIntervalAffine(fnc.body.get, inputValMap, inputErrorMap, uniformPrecision) + uniformRoundoff_IA_AA(fnc.body.get, inputValMap, inputErrorMap, uniformPrecision, trackRoundoffErrs) case ("affine", "affine") => - errorAffineAffine(fnc.body.get, inputValMap, inputErrorMap, uniformPrecision) + uniformRoundoff_AA_AA(fnc.body.get, inputValMap, inputErrorMap, uniformPrecision, trackRoundoffErrs) case ("smt", "affine") => - errorSMTAffine(fnc.body.get, inputValMap, inputErrorMap, uniformPrecision) + uniformRoundoff_SMT_AA(fnc.body.get, inputValMap, inputErrorMap, uniformPrecision, trackRoundoffErrs) // default is to use the method that attaches the info to trees. case ("subdiv", _) => - evaluateSubdiv( - fnc.body.get, + val tmp = doIntervalSubdivision( //evaluateSubdiv( + fnc.body.get, lang.TreeOps.freeVariablesOf(fnc.body.get), inputValMap, inputErrorMap, trackRoundoffErrs, uniformPrecision) + (tmp._2, tmp._1) case _ => reporter.fatalError(s"Your combination of $rangeMethod and $errorMethod" + "for computing ranges and errors is not supported.") + null } - - } + roundoffErrorMap += (fnc.id -> (errorMap.map(x => x._1 -> x._2.xhi))) + rangeResMap += (fnc.id -> rangeMap) + (fnc.id -> (resError, resRange)) + }).toMap timer.stop ctx.reporter.info(s"Finished $name") - (ctx, prg) + + (ctx.copy( + resultAbsoluteErrors = res.map(x => (x._1 -> x._2._1)), + resultRealRanges = res.map(x => (x._1 -> x._2._2)), + intermediateAbsoluteErrors = roundoffErrorMap, + intermediateRanges = rangeResMap), prg) } diff --git a/src/main/scala/daisy/analysis/RangePhase.scala b/src/main/scala/daisy/analysis/RangePhase.scala index 53e8f40d85fdca6bcc9991f2558d3199623af37e..3cdfdc8e7022efac68f5394f4690f9f00069d65f 100644 --- a/src/main/scala/daisy/analysis/RangePhase.scala +++ b/src/main/scala/daisy/analysis/RangePhase.scala @@ -70,7 +70,7 @@ object RangePhase extends DaisyPhase { case "smt" => evaluateSMT(fnc.body.get, Map.empty) case "subdiv" => - evaluateSubdiv(fnc.body.get, ctx.inputRanges(fnc.id), Map.empty) + evaluateSubdiv(fnc.body.get, ctx.specInputRanges(fnc.id), Map.empty) } diff --git a/src/main/scala/daisy/analysis/SpecsProcessingPhase.scala b/src/main/scala/daisy/analysis/SpecsProcessingPhase.scala index fb35b05605b8faec8a5254894cb36b698e70304f..43717feba44011f2859229d4d84f26d0543775fb 100644 --- a/src/main/scala/daisy/analysis/SpecsProcessingPhase.scala +++ b/src/main/scala/daisy/analysis/SpecsProcessingPhase.scala @@ -121,8 +121,8 @@ object SpecsProcessingPhase extends DaisyPhase { reporter.debug("range bounds: " + resRanges.mkString("\n")) reporter.debug("error bounds: " + resErrors.mkString("\n")) - (ctx.copy(inputRanges = allRanges, inputErrors = allErrors, - resultRangeBounds = resRanges, resultErrorBounds = resErrors), prg) + (ctx.copy(specInputRanges = allRanges, specInputErrors = allErrors, + specResultRangeBounds = resRanges, specResultErrorBounds = resErrors), prg) } diff --git a/src/main/scala/daisy/backend/CertificatePhase.scala b/src/main/scala/daisy/backend/CertificatePhase.scala index d929a9dae0d13b2d1a1abbb13ce0e6976ee8cf4c..63d134355e1dc0a21ea11a5c53af253750dc2d62 100644 --- a/src/main/scala/daisy/backend/CertificatePhase.scala +++ b/src/main/scala/daisy/backend/CertificatePhase.scala @@ -5,7 +5,7 @@ import daisy.lang.{ScalaPrinter, PrettyPrinter} import lang.Trees._ import lang.Identifiers._ import lang.NumAnnotation -import utils.Interval._ +import utils.Interval import utils.Rational import analysis.SpecsProcessingPhase @@ -20,7 +20,9 @@ object CertificatePhase extends DaisyPhase { def run(ctx: Context, prg: Program): (Context, Program) = { val reporter = ctx.reporter - val proverToUse = ctx.findOption(Main.optionValidators) + val prover = ctx.findOption(Main.optionValidators) + val errorMap = ctx.resultAbsoluteErrors + val rangeMap = ctx.resultRealRanges def writeToFile (fileContent:String, prover:String, fname:String){ import java.io.FileWriter @@ -30,9 +32,9 @@ object CertificatePhase extends DaisyPhase { "./coq/output/certificate_" + prg.id + "_" + fname + ".v" else if (prover == "hol4") - "./hol4/output/certificate_" + prg.id + "_" + fname + "Script.sml" - else - "./hol-light/output/certificate_" + prg.id + "_" + fname + ".hl" + "./hol4/output/certificate_" + prg.id + "_" + fname + "Script.sml" + else + "./hol-light/output/certificate_" + prg.id + "_" + fname + ".hl" val fstream = new FileWriter(fileLocation) val out = new BufferedWriter(fstream) out.write(fileContent) @@ -44,7 +46,10 @@ object CertificatePhase extends DaisyPhase { val thePrecondition = fnc.precondition val theBody = fnc.body - proverToUse match { + val errorMap = ctx.intermediateAbsoluteErrors.get(fnc.id).get + val rangeMap = ctx.intermediateRanges.get(fnc.id).get + + prover match { case Some(prv) => thePrecondition match { case Some (pre) => @@ -55,7 +60,7 @@ object CertificatePhase extends DaisyPhase { //generate the precondition val (thePreconditionFunction, functionName) = getPrecondFunction(pre, reporter, prv) //the analysis result function - val (analysisResultText, analysisResultName) = getAbsEnvDef(theBody.get, prv) + val (analysisResultText, analysisResultName) = getAbsEnvDef(theBody.get, errorMap, rangeMap, prv) //generate the final evaluation statement val functionCall = getComputeExpr(lastGenName,analysisResultName,functionName,prv) //val functionCall = getAllComputeExps(theBody.get, analysisResultName, functionName, prv) @@ -71,16 +76,11 @@ object CertificatePhase extends DaisyPhase { //write to the output file writeToFile(certificateText,prv,fnc.id.toString) - case None => reporter.fatalError ("No Precondition specified") + case None => reporter.fatalError ("No Precondition specified.") } - case None => reporter.fatalError ("Failure") + case None => reporter.fatalError ("No theorem prover specified.") } - } - //generate definitions for "program" - - //obtain analysis result from ctx? - //generate abstract environment for analysis result (ctx, prg) } @@ -377,57 +377,59 @@ object CertificatePhase extends DaisyPhase { } } - private def coqAbsEnv (e:Expr, cont:String) :String = + private def coqAbsEnv (e:Expr, errorMap:Map[Expr, Rational], rangeMap:Map[Expr, Interval], cont:String) :String = { //must be set TODO:Assert? + //since we already generated the AST, we can be sure to have generated the name + assert (expressionNames(e) != None) val nameE = expressionNames(e) e match { case x @ Variable(id) => - val intvE = coqInterval((x.interval.xlo,x.interval.xhi)) - val errE = x.absError.toFractionString.replace("/","#") + val intvE = coqInterval((rangeMap(x).xlo,rangeMap(x).xhi)) + val errE = errorMap(x).toFractionString.replace("/","#") "if expEqBool e " + nameE + " then (" + intvE + "," + errE + ")\n" + " else " + cont case x @ RealLiteral(r) => - val intvE = coqInterval((x.interval.xlo,x.interval.xhi)) - val errE = x.absError.toFractionString.replace("/","#") + val intvE = coqInterval((rangeMap(x).xlo,rangeMap(x).xhi)) + val errE = errorMap(x).toFractionString.replace("/","#") "if expEqBool e " + nameE + " then (" + intvE + "," + errE + ")\n" + " else " + cont case x @ Plus(lhs, rhs) => - val lFun = coqAbsEnv (lhs, cont) - val rFun = coqAbsEnv (rhs, lFun) - val intvE = coqInterval((x.interval.xlo,x.interval.xhi)) - val errE = x.absError.toFractionString.replace("/","#") + val lFun = coqAbsEnv (lhs, errorMap, rangeMap, cont) + val rFun = coqAbsEnv (rhs, errorMap, rangeMap, lFun) + val intvE = coqInterval((rangeMap(x).xlo,rangeMap(x).xhi)) + val errE = errorMap(x).toFractionString.replace("/","#") "if expEqBool e " + nameE + " then (" + intvE + "," + errE + ")\n" + " else " + rFun case x @ Minus(lhs, rhs) => - val lFun = coqAbsEnv (lhs, cont) - val rFun = coqAbsEnv (rhs, lFun) - val intvE = coqInterval((x.interval.xlo,x.interval.xhi)) - val errE = x.absError.toFractionString.replace("/","#") + val lFun = coqAbsEnv (lhs, errorMap, rangeMap, cont) + val rFun = coqAbsEnv (rhs, errorMap, rangeMap, lFun) + val intvE = coqInterval((rangeMap(x).xlo,rangeMap(x).xhi)) + val errE = errorMap(x).toFractionString.replace("/","#") "if expEqBool e " + nameE + " then (" + intvE + "," + errE + ")\n" + " else " + rFun case x @ Times(lhs, rhs) => - val lFun = coqAbsEnv (lhs, cont) - val rFun = coqAbsEnv (rhs, lFun) - val intvE = coqInterval((x.interval.xlo,x.interval.xhi)) - val errE = x.absError.toFractionString.replace("/","#") + val lFun = coqAbsEnv (lhs, errorMap, rangeMap, cont) + val rFun = coqAbsEnv (rhs, errorMap, rangeMap, lFun) + val intvE = coqInterval((rangeMap(x).xlo,rangeMap(x).xhi)) + val errE = errorMap(x).toFractionString.replace("/","#") "if expEqBool e " + nameE + " then (" + intvE + "," + errE + ")\n" + " else " + rFun case x @ Division(lhs, rhs) => - val lFun = coqAbsEnv (lhs, cont) - val rFun = coqAbsEnv (rhs, lFun) - val intvE = coqInterval((x.interval.xlo,x.interval.xhi)) - val errE = x.absError.toFractionString.replace("/","#") + val lFun = coqAbsEnv (lhs, errorMap, rangeMap, cont) + val rFun = coqAbsEnv (rhs, errorMap, rangeMap, lFun) + val intvE = coqInterval((rangeMap(x).xlo,rangeMap(x).xhi)) + val errE = errorMap(x).toFractionString.replace("/","#") "if expEqBool e " + nameE + " then (" + intvE + "," + errE + ")\n" + " else " + rFun @@ -436,57 +438,57 @@ object CertificatePhase extends DaisyPhase { } } - private def hol4AbsEnv (e:Expr, cont:String) :String = + private def hol4AbsEnv (e:Expr, errorMap:Map[Expr, Rational], rangeMap:Map[Expr, Interval], cont:String) :String = { //must be set TODO:Assert? val nameE = expressionNames(e) e match { case x @ Variable(id) => - val intvE = hol4Interval((x.interval.xlo,x.interval.xhi)) - val errE = x.absError.toFractionString + val intvE = hol4Interval((rangeMap(x).xlo,rangeMap(x).xhi)) + val errE = errorMap(x).toFractionString "if e = " + nameE + " then (" + intvE + "," + errE + ")\n" + " else " + cont case x @ RealLiteral(r) => - val intvE = hol4Interval((x.interval.xlo,x.interval.xhi)) - val errE = x.absError.toFractionString + val intvE = hol4Interval((rangeMap(x).xlo,rangeMap(x).xhi)) + val errE = errorMap(x).toFractionString "if e = " + nameE + " then (" + intvE + "," + errE + ")\n" + " else " + cont case x @ Plus(lhs, rhs) => - val lFun = hol4AbsEnv (lhs, cont) - val rFun = hol4AbsEnv (rhs, lFun) - val intvE = hol4Interval((x.interval.xlo,x.interval.xhi)) - val errE = x.absError.toFractionString + val lFun = hol4AbsEnv (lhs, errorMap, rangeMap, cont) + val rFun = hol4AbsEnv (rhs, errorMap, rangeMap, lFun) + val intvE = hol4Interval((rangeMap(x).xlo,rangeMap(x).xhi)) + val errE = errorMap(x).toFractionString "if e = " + nameE + " then (" + intvE + "," + errE + ")\n" + " else " + rFun case x @ Minus(lhs, rhs) => - val lFun = hol4AbsEnv (lhs, cont) - val rFun = hol4AbsEnv (rhs, lFun) - val intvE = hol4Interval((x.interval.xlo,x.interval.xhi)) - val errE = x.absError.toFractionString + val lFun = hol4AbsEnv (lhs, errorMap, rangeMap, cont) + val rFun = hol4AbsEnv (rhs, errorMap, rangeMap, lFun) + val intvE = hol4Interval((rangeMap(x).xlo,rangeMap(x).xhi)) + val errE = errorMap(x).toFractionString "if e = " + nameE + " then (" + intvE + "," + errE + ")\n" + " else " + rFun case x @ Times(lhs, rhs) => - val lFun = hol4AbsEnv (lhs, cont) - val rFun = hol4AbsEnv (rhs, lFun) - val intvE = hol4Interval((x.interval.xlo,x.interval.xhi)) - val errE = x.absError.toFractionString + val lFun = hol4AbsEnv (lhs, errorMap, rangeMap, cont) + val rFun = hol4AbsEnv (rhs, errorMap, rangeMap, lFun) + val intvE = hol4Interval((rangeMap(x).xlo,rangeMap(x).xhi)) + val errE = errorMap(x).toFractionString "if e = " + nameE + " then (" + intvE + "," + errE + ")\n" + " else " + rFun case x @ Division(lhs, rhs) => - val lFun = hol4AbsEnv (lhs, cont) - val rFun = hol4AbsEnv (rhs, lFun) - val intvE = hol4Interval((x.interval.xlo,x.interval.xhi)) - val errE = x.absError.toFractionString + val lFun = hol4AbsEnv (lhs, errorMap, rangeMap, cont) + val rFun = hol4AbsEnv (rhs, errorMap, rangeMap, lFun) + val intvE = hol4Interval((rangeMap(x).xlo,rangeMap(x).xhi)) + val errE = errorMap(x).toFractionString "if e = " + nameE + " then (" + intvE + "," + errE + ")\n" + " else " + rFun @@ -496,57 +498,57 @@ object CertificatePhase extends DaisyPhase { } } - private def holLightAbsEnv (e:Expr, cont:String) :String = + private def holLightAbsEnv (e:Expr, errorMap:Map[Expr, Rational], rangeMap:Map[Expr, Interval], cont:String) :String = { //must be set TODO:Assert? val nameE = expressionNames(e) e match { case x @ Variable(id) => val intvE = holLightInterval((x.interval.xlo,x.interval.xhi)) - val errE = x.absError.toFractionString.replace("(","(&") + val errE = errorMap(x).toFractionString.replace("(","(&") "if e = " + nameE + " then (" + intvE + "," + errE + ")\n" + " else " + cont case x @ RealLiteral(r) => val intvE = holLightInterval((x.interval.xlo,x.interval.xhi)) - val errE = x.absError.toFractionString.replace("(","(&") + val errE = errorMap(x).toFractionString.replace("(","(&") "if e = " + nameE + " then (" + intvE + "," + errE + ")\n" + " else " + cont case x @ Plus(lhs, rhs) => - val lFun = holLightAbsEnv (lhs, cont) - val rFun = holLightAbsEnv (rhs, lFun) + val lFun = holLightAbsEnv (lhs, errorMap, rangeMap, cont) + val rFun = holLightAbsEnv (rhs, errorMap, rangeMap, lFun) val intvE = holLightInterval((x.interval.xlo,x.interval.xhi)) - val errE = x.absError.toFractionString.replace("(","(&") + val errE = errorMap(x).toFractionString.replace("(","(&") "if e = " + nameE + " then (" + intvE + "," + errE + ")\n" + " else " + rFun case x @ Minus(lhs, rhs) => - val lFun = holLightAbsEnv (lhs, cont) - val rFun = holLightAbsEnv (rhs, lFun) + val lFun = holLightAbsEnv (lhs, errorMap, rangeMap, cont) + val rFun = holLightAbsEnv (rhs, errorMap, rangeMap, lFun) val intvE = holLightInterval((x.interval.xlo,x.interval.xhi)) - val errE = x.absError.toFractionString.replace("(","(&") + val errE = errorMap(x).toFractionString.replace("(","(&") "if e = " + nameE + " then (" + intvE + "," + errE + ")\n" + " else " + rFun case x @ Times(lhs, rhs) => - val lFun = holLightAbsEnv (lhs, cont) - val rFun = holLightAbsEnv (rhs, lFun) + val lFun = holLightAbsEnv (lhs, errorMap, rangeMap, cont) + val rFun = holLightAbsEnv (rhs, errorMap, rangeMap, lFun) val intvE = holLightInterval((x.interval.xlo,x.interval.xhi)) - val errE = x.absError.toFractionString.replace("(","(&") + val errE = errorMap(x).toFractionString.replace("(","(&") "if e = " + nameE + " then (" + intvE + "," + errE + ")\n" + " else " + rFun case x @ Division(lhs, rhs) => - val lFun = holLightAbsEnv (lhs, cont) - val rFun = holLightAbsEnv (rhs, lFun) + val lFun = holLightAbsEnv (lhs, errorMap, rangeMap, cont) + val rFun = holLightAbsEnv (rhs, errorMap, rangeMap, lFun) val intvE = holLightInterval((x.interval.xlo,x.interval.xhi)) - val errE = x.absError.toFractionString.replace("(","(&") + val errE = errorMap(x).toFractionString.replace("(","(&") "if e = " + nameE + " then (" + intvE + "," + errE + ")\n" + " else " + rFun @@ -556,15 +558,18 @@ object CertificatePhase extends DaisyPhase { } } - private def getAbsEnvDef (e:Expr, prv:String) :(String,String)= + private def getAbsEnvDef (e:Expr, errorMap:Map[Expr, Rational], rangeMap:Map[Expr, Interval], prv:String) :(String,String)= if (prv == "coq") - ("Definition absenv :analysisResult := \nfun (e:exp Q) =>\n"+ coqAbsEnv(e, "((0#1,0#1),0#1)") + ".", + ("Definition absenv :analysisResult := \nfun (e:exp Q) =>\n" + + coqAbsEnv(e, errorMap, rangeMap, "((0#1,0#1),0#1)") + ".", "absenv") else if (prv == "hol4") - ("val absenv_def = Define `absenv:analysisResult = \n\\e. \n" + hol4AbsEnv(e, "((0,1),1)") + "`;", + ("val absenv_def = Define `absenv:analysisResult = \n\\e. \n" + + hol4AbsEnv(e, errorMap, rangeMap, "((0,1),1)") + "`;", "absenv") else - ("let absenv = define `absenv:analysisResult = \n\\e. \n" + holLightAbsEnv(e, "((&0,&1),&1)") + "`;;", + ("let absenv = define `absenv:analysisResult = \n\\e. \n" + + holLightAbsEnv(e, errorMap, rangeMap, "((&0,&1),&1)") + "`;;", "absenv") private def getComputeExpr (lastGenName:String, analysisResultName:String,precondName:String, prover:String) :String= diff --git a/src/main/scala/daisy/backend/CodeGenerationPhase.scala b/src/main/scala/daisy/backend/CodeGenerationPhase.scala index 281721ba91863899044aaa67406d130decc750dc..7fce26044d64555e41cae980dfc6530c12c21a3b 100755 --- a/src/main/scala/daisy/backend/CodeGenerationPhase.scala +++ b/src/main/scala/daisy/backend/CodeGenerationPhase.scala @@ -129,7 +129,7 @@ object CodeGenerationPhase extends DaisyPhase { * @param fixed the (uniform) fixed-point precision to use // TODO: we also need to adjust the types, no? */ - private def toFixedPointCode(expr: Expr, fixed: Fixed): Expr = (expr: @unchecked) match { + def toFixedPointCode(expr: Expr, fixed: Fixed): Expr = (expr: @unchecked) match { case x @ Variable(id) => fixed match { case Fixed(16) => Variable(id.changeType(Int32Type)) case Fixed(32) => Variable(id.changeType(Int64Type)) diff --git a/src/main/scala/daisy/lang/TreeOps.scala b/src/main/scala/daisy/lang/TreeOps.scala index d0a61e578f7d8c0f4b7d8643d06f66d956a57401..73de3c8b5cccde7597347787e700f37fb87564b1 100644 --- a/src/main/scala/daisy/lang/TreeOps.scala +++ b/src/main/scala/daisy/lang/TreeOps.scala @@ -173,7 +173,7 @@ object TreeOps { } /** Returns the set of free variables in an expression */ - def variablesOf(expr: Expr): Set[Identifier] = { + def freeVariablesOf(expr: Expr): Set[Identifier] = { fold[Set[Identifier]] { case (e, subs) => val subvs = subs.flatten.toSet @@ -186,6 +186,20 @@ object TreeOps { }(expr) } + /** Returns the set of all variables occuring in an expression */ + def allVariablesOf(expr: Expr): Set[Identifier] = { + fold[Set[Identifier]] { + case (e, subs) => + val subvs = subs.flatten.toSet + e match { + case Variable(i) => subvs + i + //case Let(i, _, _) => subvs - i + //case Lambda(args, _) => subvs -- args.map(_.id) + case _ => subvs + } + }(expr) + } + /** A term is an expression which (for our purposes) does not contain propositional logic connectives. diff --git a/src/main/scala/daisy/solvers/Solver.scala b/src/main/scala/daisy/solvers/Solver.scala index b5678edb262c31cb15157373c1b742b7ec38c336..59a9c4f801b1c99447c229956716f83a5fff4e40 100644 --- a/src/main/scala/daisy/solvers/Solver.scala +++ b/src/main/scala/daisy/solvers/Solver.scala @@ -22,7 +22,7 @@ import lang.Types._ import lang.Trees._ import lang.Identifiers._ import utils.Bijection -import lang.TreeOps.variablesOf +import lang.TreeOps.freeVariablesOf import utils.Rational import Rational._ import lang.Constructors._ @@ -102,7 +102,7 @@ abstract class SMTLibSolver(val context: Context) { /* Add a constraint */ def assertConstraint(expr: Expr): Unit = { try { - variablesOf(expr).foreach(declareVariable) + freeVariablesOf(expr).foreach(declareVariable) val term = toSMT(expr)(Map()) emit(SMTAssert(term)) diff --git a/src/main/scala/daisy/utils/DynamicEvaluators.scala b/src/main/scala/daisy/utils/DynamicEvaluators.scala new file mode 100644 index 0000000000000000000000000000000000000000..4d7646b9bc8d25ddbbe1417ff3fbcb8e55810287 --- /dev/null +++ b/src/main/scala/daisy/utils/DynamicEvaluators.scala @@ -0,0 +1,149 @@ + +package daisy +package utils + +import analysis._ +import analysis.DynamicPhase._ +import lang.Trees._ +import lang.Identifiers._ +import lang.NumAnnotation +import Interval._ +import Rational._ +import FinitePrecision._ + +import daisy.analysis.Sampler._ +import MPFRFloat.{abs => mpfr_abs, max => mpfr_max, min => mpfr_min} + +import scala.collection.immutable._ + +trait DynamicEvaluators { + + var dynamicSamplesDefault: Int = 100000 + + /** + * Estimates the (roundoff) error of an expression dynamically + * by executing the finite-precision computation side-by-side with a higher- + * precision one in MPFR. + * This version does not automatically add any roundoff to the inputs, + * you have to specify input errors manually, or map all variable identifiers + * to zero if you don't want any input errors. + * Errors added to each input are randomly taken from the interval + * [- input-error, + input-error] + * + * @param expr expression whose error is to be calculated + * @param inputValMap map from variable identifier to its input interval + * @param inputErrorMap map from variable identifier to its max. input error + * @param dynamicSamples how many samples to use + */ + def errorDynamic( + expr: Expr, + inputValMap: Map[Identifier, Interval], + inputErrorMap: Map[Identifier, Rational], + dynamicSamples: Int = dynamicSamplesDefault): Rational = { + + val inputRanges: Map[Identifier, Interval] = inputValMap.map({ + case (id, i) => + (id, Interval(i.mid - inputRangeFactor * i.radius, + i.mid + inputRangeFactor * i.radius)) + }) + + val rand = new util.Random(3691285) + val sampler = new Uniform(inputRanges, 485793) //no seed = System millis + // TODO: this should not be using the ErrorMeasurer as it only needs the abs error + val measurer = new ErrorMeasurerMPFR() + var currentMaxAbsMPFR = measurer.maxAbsError + var currentMaxAbs: Double = measurer.maxAbsError.doubleValue + + var i = 0 + while (i < dynamicSamples) { + i = i + 1 + // no input errors + val dblInputsWithoutErrors: Map[Identifier, Double] = sampler.next + // baseline without errors + val mpfrInputs: Map[Identifier, MPFRFloat] = dblInputsWithoutErrors.map({ + case (x, d) => (x -> MPFRFloat.fromDouble(d)) + }) + // add errors to lower precision version, if there are any + val dblInputs: Map[Identifier, Double] = dblInputsWithoutErrors.map({ + case (x, d) => + if (rand.nextBoolean) (x -> (d + inputErrorMap(x).toDouble * rand.nextDouble)) + else (x -> (d - inputErrorMap(x).toDouble * rand.nextDouble)) + }) + + val dblOutput: Double = evalDouble(expr, dblInputs) + val mpfrOutput: MPFRFloat = evalMPFR(expr, mpfrInputs) + + measurer.nextValues(dblOutput, mpfrOutput) + + // Invariant that absolute errors have to grow monotonically + assert(currentMaxAbsMPFR <= measurer.maxAbsError) + currentMaxAbsMPFR = measurer.maxAbsError + + assert(currentMaxAbs <= measurer.maxAbsError.doubleValue) + currentMaxAbs = measurer.maxAbsError.doubleValue + + } + Rational.fromString(measurer.maxAbsError.toString()) + + } + + + + /** + * Estimates the (roundoff) error of an expression dynamically + * by executing the finite-precision computation side-by-side with a higher- + * precision one in MPFR. + * In this version roundoff errors on inputs are considered by converting + * a double-valued sample into a String and considering that string as the + * 'baseline'. If you don't want roundoff errors, use the other errorDynamic + * function. + * + * @param expr expression whose error is to be calculated + * @param inputValMap map from variable identifier to its input interval + * @param dynamicSamples how many samples to use + */ + def errorDynamicWithInputRoundoff( + expr: Expr, + inputValMap: Map[Identifier, Interval], + dynamicSamples: Int = dynamicSamplesDefault): Rational = { + + val inputRanges: Map[Identifier, Interval] = inputValMap.map({ + case (id, i) => + (id, Interval(i.mid - inputRangeFactor * i.radius, + i.mid + inputRangeFactor * i.radius)) + }) + + val sampler = new Uniform(inputRanges, 485793) //no seed = System millis + // TODO: no need for the ErrorMeasurer + val measurer = new ErrorMeasurerMPFR() + var currentMaxAbsMPFR = measurer.maxAbsError + var currentMaxAbs: Double = measurer.maxAbsError.doubleValue + + var i = 0 + while (i < dynamicSamples) { + i = i + 1 + // roundoff errors on inputs + val dblInputs: Map[Identifier, Double] = sampler.next + val mpfrInputs: Map[Identifier, MPFRFloat] = dblInputs.map({ + case (x, d) => (x -> MPFRFloat.fromString(d.toString)) + }) + + val dblOutput: Double = evalDouble(expr, dblInputs) + val mpfrOutput: MPFRFloat = evalMPFR(expr, mpfrInputs) + + measurer.nextValues(dblOutput, mpfrOutput) + + // Invariant that absolute errors have to grow monotonically + assert(currentMaxAbsMPFR <= measurer.maxAbsError) + currentMaxAbsMPFR = measurer.maxAbsError + + assert(currentMaxAbs <= measurer.maxAbsError.doubleValue) + currentMaxAbs = measurer.maxAbsError.doubleValue + + } + Rational.fromString(measurer.maxAbsError.toString()) + + } + + +} \ No newline at end of file diff --git a/src/main/scala/daisy/utils/ErrorFunctions.scala b/src/main/scala/daisy/utils/ErrorFunctions.scala deleted file mode 100644 index 3ac9d48aa961cc290294a5792681b1c76c5bbbd9..0000000000000000000000000000000000000000 --- a/src/main/scala/daisy/utils/ErrorFunctions.scala +++ /dev/null @@ -1,867 +0,0 @@ - -package daisy -package utils - -import analysis._ -import analysis.DynamicPhase._ -import lang.Trees._ -import lang.Identifiers._ -import lang.NumAnnotation -import Interval._ -import Rational._ -import FinitePrecision._ - -import daisy.analysis.Sampler._ -import MPFRFloat.{abs => mpfr_abs, max => mpfr_max, min => mpfr_min} - -import scala.collection.immutable._ - -/** - * Implements functions for calculating error - * Given an expression, its input value map, - * and its initial errors - * Parameters can be tweaked by the implementing object - */ -trait ErrorFunctions { - - var trackRoundoffErrs: Boolean = true - var attachToTree: Boolean = false - var dynamicSamplesDefault: Int = 100000 - //FIXME: Why does this not compile without this? - //var defaultUniformPrecision: Precision = Float64 - var uniformPrecision: Precision = Float64 - - var reporter: Reporter - implicit val debugSection: DebugSection - - - /** - * Calculate static error using interval ranges and interval errors - * - * @param expr Expression whose erros is to be calculated - * @param inputValMap Map from function identifier to its input value interval - * @param inputErrorMap Map from function - */ - def errorIntervalInterval( - expr: Expr, - inputValMap: Map[Identifier, Interval], - inputErrorMap: Map[Identifier, Rational]): Rational = { - val (_, error) = evaluate[Interval, Interval]( - expr, - inputValMap, - inputErrorMap.map(x => (x._1 -> Interval.fromError(x._2))), - Interval.apply, - Interval.zero, - Interval.fromError, - Interval.apply, - attachToTree, - trackRoundoffErrs, - uniformPrecision - ) - maxAbs(error.toInterval) - } - - /** - * Calculates static error using interval ranges - * and affine errors - * - * @param expr Expression whose error is to be calculated - * @param inputValMap Map from function identifier to its input value interval - * @param inputErrorMap Map from function - */ - def errorIntervalAffine( - expr: Expr, - inputValMap: Map[Identifier, Interval], - inputErrorMap: Map[Identifier, Rational], - uniformPrecision: Precision): Rational = { - val (_, error) = evaluate[Interval, AffineForm]( - expr, - inputValMap, - inputErrorMap.map(x => (x._1 -> AffineForm.fromError(x._2))), - Interval.apply, - AffineForm.zero, - AffineForm.fromError, - AffineForm.apply, - attachToTree, - trackRoundoffErrs, - uniformPrecision - ) - maxAbs(error.toInterval) - } - - /** - * Calculates static error using affine ranges - * and affine errors - * - * @param expr Expression whose error is to be calculated - * @param inputValMap Map from function identifier to its input value interval - * @param inputErrorMap Map from function - */ - def errorAffineAffine( - expr: Expr, - inputValMap: Map[Identifier, Interval], - inputErrorMap: Map[Identifier, Rational], - uniformPrecision: Precision): Rational = { - val (_, error) = evaluate[AffineForm, AffineForm]( - expr, - inputValMap.map(x => (x._1 -> AffineForm(x._2))), - inputErrorMap.map(x => (x._1 -> AffineForm.fromError(x._2))), - AffineForm.apply, - AffineForm.zero, - AffineForm.fromError, - (x: AffineForm) => x, // Affine => Affine - attachToTree, - trackRoundoffErrs, - uniformPrecision - ) - maxAbs(error.toInterval) - } - - /** - * Calculates static error using SMT ranges - * and affine errors - * - * @param expr Expression whose error is to be calculated - * @param inputValMap Map from function identifier to its input value interval - * @param inputErrorMap Map from function - */ - def errorSMTAffine( - expr: Expr, - inputValMap: Map[Identifier, Interval], - inputErrorMap: Map[Identifier, Rational], - uniformPrecision: Precision): Rational = { - val (_, error) = evaluate[SMTRange, AffineForm]( - expr, - inputValMap.map({ case (id, int) => (id -> SMTRange(Variable(id), int)) }), - inputErrorMap.map(x => (x._1 -> AffineForm.fromError(x._2))), - SMTRange.apply, - AffineForm.zero, - AffineForm.fromError, - (x: SMTRange) => AffineForm(x.toInterval), // SMTRange => AffineForm - attachToTree, - trackRoundoffErrs, - uniformPrecision - ) - maxAbs(error.toInterval) - } - - /** - * Estimates the (roundoff) error of an expression dynamically - * by executing the finite-precision computation side-by-side with a higher- - * precision one in MPFR. - * This version does not automatically add any roundoff to the inputs, - * you have to specify input errors manually, or map all variable identifiers - * to zero if you don't want any input errors. - * Errors added to each input are randomly taken from the interval - * [- input-error, + input-error] - * - * @param expr expression whose error is to be calculated - * @param inputValMap map from variable identifier to its input interval - * @param inputErrorMap map from variable identifier to its max. input error - * @param dynamicSamples how many samples to use - */ - def errorDynamic( - expr: Expr, - inputValMap: Map[Identifier, Interval], - inputErrorMap: Map[Identifier, Rational], - dynamicSamples: Int = dynamicSamplesDefault): Rational = { - - val inputRanges: Map[Identifier, Interval] = inputValMap.map({ - case (id, i) => - (id, Interval(i.mid - inputRangeFactor * i.radius, - i.mid + inputRangeFactor * i.radius)) - }) - - val rand = new util.Random(3691285) - val sampler = new Uniform(inputRanges, 485793) //no seed = System millis - // TODO: this should not be using the ErrorMeasurer as it only needs the abs error - val measurer = new ErrorMeasurerMPFR() - var currentMaxAbsMPFR = measurer.maxAbsError - var currentMaxAbs: Double = measurer.maxAbsError.doubleValue - - var i = 0 - while (i < dynamicSamples) { - i = i + 1 - // no input errors - val dblInputsWithoutErrors: Map[Identifier, Double] = sampler.next - // baseline without errors - val mpfrInputs: Map[Identifier, MPFRFloat] = dblInputsWithoutErrors.map({ - case (x, d) => (x -> MPFRFloat.fromDouble(d)) - }) - // add errors to lower precision version, if there are any - val dblInputs: Map[Identifier, Double] = dblInputsWithoutErrors.map({ - case (x, d) => - if (rand.nextBoolean) (x -> (d + inputErrorMap(x).toDouble * rand.nextDouble)) - else (x -> (d - inputErrorMap(x).toDouble * rand.nextDouble)) - }) - - val dblOutput: Double = evalDouble(expr, dblInputs) - val mpfrOutput: MPFRFloat = evalMPFR(expr, mpfrInputs) - - measurer.nextValues(dblOutput, mpfrOutput) - - // Invariant that absolute errors have to grow monotonically - assert(currentMaxAbsMPFR <= measurer.maxAbsError) - currentMaxAbsMPFR = measurer.maxAbsError - - assert(currentMaxAbs <= measurer.maxAbsError.doubleValue) - currentMaxAbs = measurer.maxAbsError.doubleValue - - } - Rational.fromString(measurer.maxAbsError.toString()) - - } - - /** - * Estimates the (roundoff) error of an expression dynamically - * by executing the finite-precision computation side-by-side with a higher- - * precision one in MPFR. - * In this version roundoff errors on inputs are considered by converting - * a double-valued sample into a String and considering that string as the - * 'baseline'. If you don't want roundoff errors, use the other errorDynamic - * function. - * - * @param expr expression whose error is to be calculated - * @param inputValMap map from variable identifier to its input interval - * @param dynamicSamples how many samples to use - */ - def errorDynamicWithInputRoundoff( - expr: Expr, - inputValMap: Map[Identifier, Interval], - dynamicSamples: Int = dynamicSamplesDefault): Rational = { - - val inputRanges: Map[Identifier, Interval] = inputValMap.map({ - case (id, i) => - (id, Interval(i.mid - inputRangeFactor * i.radius, - i.mid + inputRangeFactor * i.radius)) - }) - - val sampler = new Uniform(inputRanges, 485793) //no seed = System millis - // TODO: no need for the ErrorMeasurer - val measurer = new ErrorMeasurerMPFR() - var currentMaxAbsMPFR = measurer.maxAbsError - var currentMaxAbs: Double = measurer.maxAbsError.doubleValue - - var i = 0 - while (i < dynamicSamples) { - i = i + 1 - // roundoff errors on inputs - val dblInputs: Map[Identifier, Double] = sampler.next - val mpfrInputs: Map[Identifier, MPFRFloat] = dblInputs.map({ - case (x, d) => (x -> MPFRFloat.fromString(d.toString)) - }) - - val dblOutput: Double = evalDouble(expr, dblInputs) - val mpfrOutput: MPFRFloat = evalMPFR(expr, mpfrInputs) - - measurer.nextValues(dblOutput, mpfrOutput) - - // Invariant that absolute errors have to grow monotonically - assert(currentMaxAbsMPFR <= measurer.maxAbsError) - currentMaxAbsMPFR = measurer.maxAbsError - - assert(currentMaxAbs <= measurer.maxAbsError.doubleValue) - currentMaxAbs = measurer.maxAbsError.doubleValue - - } - Rational.fromString(measurer.maxAbsError.toString()) - - } - - /** - * Analyses ranges and errors with different methods. - * - * To check which combinations are definitely supported, see the run method. - * There you can also see the default values for the (good few) params needed. - * - * @param expr expression to be evaluated - * @param initValMap initial ranges for all free variables of expr - * @param initErrorMap errors (incl. roundoff errors) for all free variables - * @param rangeFromReal constructor to create a range from a rational with - * the range arithmetic chosed to track ranges - * @param zeroError expression denoting no error in the range arithmetic - * chose for tracking errors - * @param fromError constructor to create the range arithm. expression from - * an (roundoff) error - * @param t2s transformation function to cast the range arith. used for tracking - * ranges into the one tracking errors - * @param attachToTree if true, ranges and errors will be attached to the tree - * Note: this can only be done once, and this method does not check - * whether a range is already attched (will crash in that case) - * @param trackRoundoffErrs if true, roundoff errors will be tracked, otherwise - * only the initial errors will be propagated - * @param uniformPrecision the default precision to use when doing absolute - * errors calculations - */ - def evaluate[T <: RangeArithmetic[T], S <: RangeArithmetic[S]]( - expr: Expr, - initValMap: Map[Identifier, T], - initErrorMap: Map[Identifier, S], - rangeFromReal: Rational => T, - zeroError: S, - fromError: Rational => S, - t2s: T => S, - attachToTree: Boolean, - trackRoundoffErrs: Boolean, - uniformPrecision: Precision): (T, S) = { - - @inline - def computeNewError(range: T, propagatedError: S): S = - if (trackRoundoffErrs) { - val actualRange: Interval = range.toInterval + propagatedError.toInterval - val rndoff = uniformPrecision.absRoundoff(actualRange) - propagatedError +/- rndoff - } else { - propagatedError - } - - def eval(e: Expr, valMap: Map[Identifier, T], errorMap: Map[Identifier, S]): (T, S) = - (e: @unchecked) match { - - case x @ Variable(id) => - val range = valMap(id) - val error = errorMap(id) - if (attachToTree) { - x.interval = range.toInterval - x.absError = maxAbs(error.toInterval) - } - //print state - reporter.debug(s"$x: ${x.interval},${x.absError}") - (range, error) - - case x @ RealLiteral(r) => - val range = rangeFromReal(r) - // was: isExactInFloats(r, uniformPrecision) || !trackRoundoffErrs) - val error = if (!trackRoundoffErrs) { - zeroError - } else { - fromError(uniformPrecision.absRoundoff(r)) - } - if (attachToTree) { - x.absError = maxAbs(error.toInterval) - x.interval = range.toInterval - } - //print state - reporter.debug(s"$x: ${x.interval},${x.absError}") - (range, error) - - case x @ Plus(lhs, rhs) => - val (rangeLHS, errorLHS) = eval(lhs, valMap, errorMap) - val (rangeRHS, errorRHS) = eval(rhs, valMap, errorMap) - - val range = rangeLHS + rangeRHS - - val propagatedError = errorLHS + errorRHS - val newError = computeNewError(range, propagatedError) - - if (attachToTree) { - x.absError = maxAbs(newError.toInterval) - x.interval = range.toInterval - } - //print state - reporter.debug(s"$x: ${x.interval},${x.absError}") - (range, newError) - - case x @ Minus(lhs, rhs) => - val (rangeLHS, errorLHS) = eval(lhs, valMap, errorMap) - val (rangeRHS, errorRHS) = eval(rhs, valMap, errorMap) - - val range = rangeLHS - rangeRHS - - val propagatedError = errorLHS - errorRHS - val newError = computeNewError(range, propagatedError) - - if (attachToTree) { - x.absError = maxAbs(newError.toInterval) - x.interval = range.toInterval - } - //print state - reporter.debug(s"$x: ${x.interval},${x.absError}") - (range, newError) - - case x @ Times(lhs, rhs) => - val (rangeLHS, errorLHS) = eval(lhs, valMap, errorMap) - val (rangeRHS, errorRHS) = eval(rhs, valMap, errorMap) - - // new range - val range: T = rangeLHS * rangeRHS - - // error propagation - val leftRangeInS = t2s(rangeLHS) - val rightRangeInS = t2s(rangeRHS) - val propagatedError = leftRangeInS * errorRHS + rightRangeInS * errorLHS + - errorLHS*errorRHS - // new roundoff - val newError = computeNewError(range, propagatedError) - - if (attachToTree) { - x.absError = maxAbs(newError.toInterval) - x.interval = range.toInterval - } - //print state - reporter.debug(s"$x: ${x.interval},${x.absError}") - (range, newError) - - case x @ Division(lhs, rhs) => - val (rangeLHS, errorLHS) = eval(lhs, valMap, errorMap) - val (rangeRHS, errorRHS) = eval(rhs, valMap, errorMap) - - // new range - val range = rangeLHS / rangeRHS - - // inverse, i.e. we are computing x * (1/y) - val rightInterval = rangeRHS.toInterval + errorRHS.toInterval // the actual interval, incl errors - - //the actual error interval can now contain 0. - //check this - if (rightInterval.xlo <= 0 && rightInterval.xhi >= 0) - throw DivisionByZeroException("trying to divide by error interval containing 0") - - val a = min(abs(rightInterval.xlo), abs(rightInterval.xhi)) - //one = 1 - val errorMultiplier: Rational = -one / (a*a) - val invErr = errorRHS * errorMultiplier - reporter.debug(s"${errorRHS}, $invErr") - // error propagation - val leftRangeInS = t2s(rangeLHS) - val inverse: T = rangeRHS.inverse - val kRangeInS = t2s(inverse) - var propagatedError = leftRangeInS * invErr + kRangeInS * errorLHS + errorLHS * invErr - reporter.debug(s"$propagatedError") - // new roundoff - val newError = computeNewError(range, propagatedError) - - if (attachToTree) { - x.absError = maxAbs(newError.toInterval) - x.interval = range.toInterval - } - //print state - reporter.debug(s"$x: ${x.interval},${x.absError}") - (range, newError) - - case x @ UMinus(t) => - val (range, error) = eval(t, valMap, errorMap) - // negation does not affect error magnitude - (-range, -error) - - case x @ Sqrt(t) => - // TODO: Not supported for fixed-points, add exception - //case FPPrecision(_) => throw UnsupportedRealFragmentException("Sqrt not supported for fixed-points.") - val (rangeT, errorT) = eval(t, valMap, errorMap) - val range = rangeT.squareRoot - - // propagated existing errors - val tInterval = rangeT.toInterval - val a = min(abs(tInterval.xlo), abs(tInterval.xhi)) - val errorMultiplier = Rational(1L, 2L) / sqrtDown(a) - val propagatedError = errorT * errorMultiplier - - // new roundoff - val newError = computeNewError(range, propagatedError) - - if (attachToTree) { - x.absError = maxAbs(newError.toInterval) - x.interval = range.toInterval - } - (range, newError) - - case Let(id, value, body) => - val (valueRange, valueError) = eval(value, valMap, errorMap) - val newMapRange = valMap + (id -> valueRange) - val newMapError = errorMap + (id -> valueError) - eval(body, newMapRange, newMapError) - } - eval(expr, initValMap, initErrorMap) - } - - - /** - * Performs an interval subdivision, recording only the result's range - * and error. - * - * The tree does not get annotated. - * The initErrorMap only contains initial errors, NOT roundoffs. - * - * @param expr expression to be analysed - * @param inputParams free variables of the expr - * @param initValMap ranges of free variables - * @param initErrorMap initial errors of free variables (if any) - * These do NOT include roundoff errors - * @param trackRoundoffErrs whether to track roundoff errors or to only - * propagate initial errors (as provided) - * @param uniformPrecision the default precision to use when doing absolute - * errors calculations - */ - def doIntervalSubdivision( - expr: Expr, - inputParams: Set[Identifier], - initValMap: Map[Identifier, Interval], - initErrorMap: Map[Identifier, Rational], - //trackInitialErrs: Boolean, - trackRoundoffErrs: Boolean, - uniformPrecision: Precision): (Interval, Rational) = { - - // TODO: this should not be hardcoded - val numSplits = 10 - - val inputsSubdiv: Seq[Map[Identifier, Interval]] = initValMap.foldLeft(Seq(Map[Identifier, Interval]()))({ - case (currSeq: Seq[Map[Identifier, Interval]], (id, intrvl)) => - val xlo = intrvl.xlo - val splitWidth = intrvl.width / numSplits - val splits: Seq[Interval] = (0 until numSplits).map(i => - Interval(xlo + i * splitWidth, xlo + (i+1) * splitWidth)) - - currSeq.flatMap( m => - splits.map(i => m + (id -> i)) - ) - }) - - val allIDs = inputParams - val missingIDs = allIDs -- initErrorMap.keySet - - // the initErrorMap will simply be empty, if initial errors are not being tracked - def getInputErrors(rangeMap: Map[Identifier, Interval]): Map[Identifier, Rational] = - if (trackRoundoffErrs){ - - initErrorMap ++ missingIDs.map(id => (id -> uniformPrecision.absRoundoff(rangeMap(id)))) - - } else { - - initErrorMap ++ missingIDs.map( id => (id -> Rational.zero)) - - } - - val (resInterval, resError) = evaluate[Interval, AffineForm]( - expr, inputsSubdiv.head, - getInputErrors(inputsSubdiv.head).map(x => (x._1 -> AffineForm.fromError(x._2))), - Interval.apply, AffineForm.zero, AffineForm.fromError, AffineForm.apply, - false, trackRoundoffErrs, uniformPrecision) - - var maxInterval = resInterval - var maxError = maxAbs(resError.toInterval) - - for (m <- inputsSubdiv.tail) { - - val (tmpRange, tmpError) = evaluate[Interval, AffineForm]( - expr, m, - getInputErrors(m).map(x => (x._1 -> AffineForm.fromError(x._2))), - Interval.apply, AffineForm.zero, AffineForm.fromError, AffineForm.apply, - false, trackRoundoffErrs, uniformPrecision) - - maxInterval = maxInterval.union(tmpRange) - maxError = max(maxError, maxAbs(tmpError.toInterval)) - } - (maxInterval, maxError) - } - - - /** - * Uses interval subdivision to compute ranges and affine arithmetic for errors. - * - * Unlike the other subdivision method, this one annotates the trees with - * information. For this, ranges are computed first with interval subdivision, - * and errors are computed in a second step. Due to this separation, the computed - * errors may not be the tightest possible. - * This methods is suitable for fixed-point arithmetic. - * - * @param expr expression to be analysed - * @param inputRanges initial ranges of free variables - * @param varErrorMap input errors (incl. roundoff) of free variables - * @param trackRoundoffErrs whether or not to track roundoff errors or to - * propagate only input errors (as given in varErrorMap) - * @param uniformPrecision the default precision to use when doing absolute - * errors calculations - */ - def evaluateSubdiv( - expr: Expr, - inputRanges: Map[Identifier, Interval], - varErrorMap: Map[Identifier, Rational], - trackRoundoffErrs: Boolean, - uniformPrecision: Precision): (Interval, Rational) = { - - @inline - def computeNewError(range: Interval, propagatedError: AffineForm): AffineForm = - if (trackRoundoffErrs) { - val actualRange: Interval = range + propagatedError.toInterval - val rndoff = uniformPrecision.absRoundoff(actualRange) - propagatedError +/- rndoff - } else { - propagatedError - } - - val numSplits = 10 - - //val inputRanges: Map[Identifier, Interval] = context.inputRanges(fncId) - - val inputsSubdiv: Seq[Map[Identifier, Interval]] = inputRanges.foldLeft(Seq(Map[Identifier, Interval]()))({ - case (currSeq: Seq[Map[Identifier, Interval]], (id, intrvl)) => - val xlo = intrvl.xlo - val splitWidth = intrvl.width / numSplits - val splits: Seq[Interval] = (0 until numSplits).map(i => - Interval(xlo + i * splitWidth, xlo + (i+1) * splitWidth)) - - currSeq.flatMap( m => - splits.map(i => m + (id -> i)) - ) - }) - - // map of maximum seen ranges of all intermediate expressions - var currentRanges: Map[Expr, Interval] = Map() - - // there is probably a helper function for this - def attachInterval(e: Expr): Unit = (e: @unchecked) match { - - case x @ Variable(id) => - if (!x.hasInterval) { - x.interval = currentRanges(x) - } - return - - case x @ RealLiteral(r) => - x.interval = Interval(r) - return - - case x @ Plus(lhs, rhs) => - attachInterval(lhs) - attachInterval(rhs) - x.interval = currentRanges(x) - - case x @ Minus(lhs, rhs) => - attachInterval(lhs) - attachInterval(rhs) - x.interval = currentRanges(x) - - case x @ Times(lhs, rhs) => - attachInterval(lhs) - attachInterval(rhs) - x.interval = currentRanges(x) - - case x @ Division(lhs, rhs) => - attachInterval(lhs) - attachInterval(rhs) - x.interval = currentRanges(x) - - case x @ UMinus(t) => - attachInterval(t) - x.interval = currentRanges(x) - - case x @ Sqrt(t) => - attachInterval(t) - x.interval = currentRanges(x) - - // TODO: we may be missing intervals on the let-defined vars - case x @ Let(id, value, body) => - attachInterval(value) - attachInterval(body) - // the range of the let is the range of the final result returned - x.interval = currentRanges(body) - currentRanges = currentRanges + (x -> x.interval) - } - - // assumes that the ranges are attached - def computeAndAttachError(e: Expr, errorMap: Map[Identifier, AffineForm]): AffineForm = (e: @unchecked) match { - case x @ Variable(id) => - val err = errorMap(id) - x.absError = maxAbs(err.toInterval) - err - - case x @ RealLiteral(r) => - val err = if (isExactInFloats(r, uniformPrecision)) { - AffineForm.zero - } else { - AffineForm.fromError(uniformPrecision.absRoundoff(r)) - } - x.absError = maxAbs(err.toInterval) - err - - case x @ Plus(lhs, rhs) => - val lErr = computeAndAttachError(lhs, errorMap) - val rErr = computeAndAttachError(rhs, errorMap) - - val propagatedError = lErr + rErr - val newError = computeNewError(x.interval, propagatedError) - - x.absError = maxAbs(newError.toInterval) - newError - - case x @ Minus(lhs, rhs) => - val lErr = computeAndAttachError(lhs, errorMap) - val rErr = computeAndAttachError(rhs, errorMap) - - val propagatedError = lErr + rErr - val newError = computeNewError(x.interval, propagatedError) - - x.absError = maxAbs(newError.toInterval) - newError - - case x @ Times(lhs, rhs) => - // propagated error - val lAA = AffineForm(lhs.asInstanceOf[NumAnnotation].interval) - val rAA = AffineForm(rhs.asInstanceOf[NumAnnotation].interval) - val lErr = computeAndAttachError(lhs, errorMap) - val rErr = computeAndAttachError(rhs, errorMap) - - var propagatedError = lAA*rErr + rAA*lErr + lErr*rErr - - // new error - val newError = computeNewError(x.interval, propagatedError) - - x.absError = maxAbs(newError.toInterval) - newError - - case x @ Division(lhs, rhs) => - // propagated error - - val rhsInterval = rhs.asInstanceOf[NumAnnotation].interval - val rErr = computeAndAttachError(rhs, errorMap) //yErr - - val rInt = rhsInterval + rErr.toInterval // the actual interval, incl errors - val a = min(abs(rInt.xlo), abs(rInt.xhi)) - val errorMultiplier = -one / (a*a) - - val invErr = rErr * AffineForm(errorMultiplier) - val lAA = AffineForm(lhs.asInstanceOf[NumAnnotation].interval) - val inverse: Interval = rhsInterval.inverse - val kAA = AffineForm(inverse) - - val lErr = computeAndAttachError(lhs, errorMap) //xErr - - // multiplication with inverse - var propagatedError = lAA*invErr + kAA*lErr + lErr*invErr - - // new error - val newError = computeNewError(x.interval, propagatedError) - - x.absError = maxAbs(newError.toInterval) - newError - - case x @ UMinus(t) => - val newError = - computeAndAttachError(t, errorMap) // no additional roundoff error - x.absError = maxAbs(newError.toInterval) - newError - - case x @ Sqrt(t) => - - // TODO: Not supported for fixed-points, add exception - //case FPPrecision(_) => throw UnsupportedRealFragmentException("Sqrt not supported for fixed-points.") - - val tInterval = t.asInstanceOf[NumAnnotation].interval - val tError = computeAndAttachError(t, errorMap) - - // propagated existing errors - val a = min(abs(tInterval.xlo), abs(tInterval.xhi)) - val errorMultiplier = Rational(1L, 2L) / sqrtDown(a) - val propagatedError = tError * AffineForm(errorMultiplier) - - // new roundoff - val newError = computeNewError(x.interval, propagatedError) - - x.absError = maxAbs(newError.toInterval) - newError - - - case Let(id, value, body) => - val aform = computeAndAttachError(value, errorMap) - val newErrorMap = errorMap + (id -> aform) - computeAndAttachError(body, newErrorMap) - } - - // collect possible intervals - // first intervals need to populate the map - evalInterval(expr, inputsSubdiv.head, - (x: Expr, i: Interval) => currentRanges = currentRanges + (x -> i)) - - // the next intervals update the map - inputsSubdiv.tail.foreach( subdiv => - evalInterval(expr, subdiv, - (x: Expr, i: Interval) => currentRanges = currentRanges + (x -> currentRanges(x).union(i)))) - - currentRanges = currentRanges ++ inputRanges.map(x => (Variable(x._1) -> x._2)) - - // traverse the tree once more and attach the correct ranges - attachInterval(expr) - - // Now that we have the ranges, we can compute the errors - val resError = computeAndAttachError(expr, varErrorMap.map(x => (x._1 -> AffineForm.fromError(x._2)))) - - (currentRanges(expr), maxAbs(resError.toInterval)) - } - - /** - * Generic interval evaluation function, which can additionally perform - * some action at every node. - * - * The action can, for example, be used to record the ranges for all - * intermediate nodes. Alternatively, it can also be a no-op. - * - * @param e expression to be analysed - * @param valMap ranges of free variables - * @param action action to be performed at each node, the node and its range - * is passed to this function - */ - def evalInterval( - e: Expr, - valMap: Map[Identifier, Interval], - action: (Expr, Interval) => Unit): Interval = (e: @unchecked) match { - - case x @ Variable(id) => - valMap(id) - - case x @ RealLiteral(r) => - Interval(r) - - case x @ Plus(lhs, rhs) => - val intrvl = evalInterval(lhs, valMap, action) + evalInterval(rhs, valMap, action) - action(x, intrvl) - //currentRanges = currentRanges + (x -> currentRanges(x).union(intrvl)) - intrvl - - case x @ Minus(lhs, rhs) => - val intrvl = evalInterval(lhs, valMap, action) - evalInterval(rhs, valMap, action) - action(x, intrvl) - //currentRanges = currentRanges + (x -> currentRanges(x).union(intrvl)) - intrvl - - case x @ Times(lhs, rhs) => - val intrvl = evalInterval(lhs, valMap, action) * evalInterval(rhs, valMap, action) - action(x, intrvl) - //currentRanges = currentRanges + (x -> currentRanges(x).union(intrvl)) - intrvl - - case x @ Division(lhs, rhs) => - try { - val intrvl = evalInterval(lhs, valMap, action) / evalInterval(rhs, valMap, action) - action(x, intrvl) - //currentRanges = currentRanges + (x -> currentRanges(x).union(intrvl)) - intrvl - } catch { - case e: utils.DivisionByZeroException => - reporter.fatalError(x.getPos, "possible division by zero") - } - - case x @ UMinus(t) => - val intrvl = - evalInterval(t, valMap, action) - action(x, intrvl) - //currentRanges = currentRanges + (x -> currentRanges(x).union(intrvl)) - intrvl - - case x @ Sqrt(t) => - try { - val intrvl = evalInterval(t, valMap, action).squareRoot - action(x, intrvl) - //currentRanges = currentRanges + (x -> currentRanges(x).union(intrvl)) - intrvl - } catch { - case e: utils.NegativeSqrtException => - reporter.fatalError(x.getPos, "possible negative square root") - } - - case x @ Let(id, value, body) => - val intrvl = evalInterval(value, valMap, action) - action(Variable(id), intrvl) - //currentRanges = currentRanges + (Variable(id) -> currentRanges(Variable(id)).union(intrvl)) - val newValMap = valMap + (id -> intrvl) - evalInterval(body, newValMap, action) - } - -} diff --git a/src/main/scala/daisy/utils/FinitePrecision.scala b/src/main/scala/daisy/utils/FinitePrecision.scala index 0b73fe0ea75dd4ad251bb9f0d72b9c40a0f4a6d1..15947441eeee8372bdb217a09d84669b25b029ae 100644 --- a/src/main/scala/daisy/utils/FinitePrecision.scala +++ b/src/main/scala/daisy/utils/FinitePrecision.scala @@ -15,41 +15,57 @@ object FinitePrecision { // if it's an integer, it's definitely representable if (r.isWhole) { true - } else if (prec == DoubleDouble || prec == QuadDouble) { - // don't want to deal with those huge number here right now, - // so return the safe answer - false } else { - val nominator = r.n.abs - val denominator = r.d.abs - - val (nomBound, denomBound) = (prec: @unchecked) match { - case Float32 => - // 2^23 - 1, 2^8 -1 - (8388607l, 255l) - case Float64 => - // 2^52 - 1, 2^11 -1 - (4503599627370496l, 2047l) - } + prec match { + case DoubleDouble | QuadDouble => false + case Fixed(_) => false + case Float32 | Float64 => + + val nominator = r.n.abs + val denominator = r.d.abs + + val (nomBound, denomBound) = (prec: @unchecked) match { + case Float32 => + // 2^23 - 1, 2^8 -1 + (8388607l, 255l) + case Float64 => + // 2^52 - 1, 2^11 -1 + (4503599627370496l, 2047l) + } - if (nominator <= nomBound && denominator <= denomBound) { - val exponent: Double = math.log(denominator.toDouble) / math.log(2) + if (nominator <= nomBound && denominator <= denomBound) { + val exponent: Double = math.log(denominator.toDouble) / math.log(2) + + if (exponent.isWhole) { + // maybe the log computations isn't sound due to roundoffs, let's sanity check: + assert(math.pow(2, exponent) == denominator) + true + } else { + false + } - if (exponent.isWhole) { - // maybe the log computations isn't sound due to roundoffs, let's sanity check: - assert(math.pow(2, exponent) == denominator) - true } else { false } + } + } + } + private def allPrec: List[Precision] = List(Float32, Float64, DoubleDouble, QuadDouble) + + def getUpperBound(lhs: Precision, rhs: Precision): Precision = (lhs, rhs) match { + case (Fixed(a), Fixed(b)) if (a == b) => lhs + case (Fixed(a), Fixed(b)) => + throw new Exception("mixed-precision currently unsupported for fixed-points") + case _ => + if (allPrec.indexOf(lhs) <= allPrec.indexOf(rhs)) { + rhs } else { - false + lhs } - } } - sealed abstract class Precision { + sealed abstract class Precision extends Ordered[Precision] { /* The range of values that are representable by this precision. */ def range: (Rational, Rational) @@ -80,6 +96,12 @@ object FinitePrecision { Rational.fromDouble(math.ulp(1.0)/2)*Rational.abs(r) //Rational.fromDouble(math.ulp(Rational.abs(r).floatValue)/2) } + + def compare(that: Precision): Int = that match { + case Float32 => 0 + case Float64 => -1 + case DoubleDouble => -1 + } } case object Float64 extends Precision { @@ -97,6 +119,12 @@ object FinitePrecision { Rational.fromDouble(math.ulp(1.0)/2)*Rational.abs(r) //Rational.fromDouble(math.ulp(Rational.abs(r).doubleValue)/2) } + + def compare(that: Precision): Int = that match { + case Float32 => 1 + case Float64 => 0 + case DoubleDouble => -1 + } } case object DoubleDouble extends Precision { @@ -115,6 +143,12 @@ object FinitePrecision { def absRoundoff(r: Rational): Rational = { doubleDoubleEps * Rational.abs(r) } + + def compare(that: Precision): Int = that match { + case Float32 => 1 + case Float64 => 1 + case DoubleDouble => 0 + } } case object QuadDouble extends Precision { @@ -133,6 +167,7 @@ object FinitePrecision { def absRoundoff(r: Rational): Rational = { quadDoubleEps * Rational.abs(r) } + def compare(that: Precision): Int = ??? } /* @@ -170,6 +205,9 @@ object FinitePrecision { // TODO: don't we have to also subtract 1 for the sign? 32 - Integer.numberOfLeadingZeros(value) } + def compare(that: Precision): Int = that match { + case Fixed(x) => bitlength.compare(x) + } } } diff --git a/src/main/scala/daisy/utils/IntervalSubdivision.scala b/src/main/scala/daisy/utils/IntervalSubdivision.scala new file mode 100644 index 0000000000000000000000000000000000000000..9e4fe92b17fe10ad110e5360a2087b4ea6ea5d71 --- /dev/null +++ b/src/main/scala/daisy/utils/IntervalSubdivision.scala @@ -0,0 +1,341 @@ + +package daisy +package utils + +import analysis._ +import analysis.DynamicPhase._ +import lang.Trees._ +import lang.Identifiers._ +import lang.NumAnnotation +import Interval._ +import Rational._ +import FinitePrecision._ + +import daisy.analysis.Sampler._ +import MPFRFloat.{abs => mpfr_abs, max => mpfr_max, min => mpfr_min} + +import scala.collection.immutable._ + +trait IntervalSubdivision extends RoundoffEvaluators { + + /** + * Performs an interval subdivision, recording only the result's range + * and error. + * + * The tree does not get annotated. + * The initErrorMap only contains initial errors, NOT roundoffs. + * + * @param expr expression to be analysed + * @param inputParams free variables of the expr + * @param initValMap ranges of free variables + * @param initErrorMap initial errors of free variables (if any) + * These do NOT include roundoff errors + * @param trackRoundoffErrs whether to track roundoff errors or to only + * propagate initial errors (as provided) + * @param uniformPrecision the default precision to use when doing absolute + * errors calculations + */ + def doIntervalSubdivision( + expr: Expr, + inputParams: Set[Identifier], + initValMap: Map[Identifier, Interval], + initErrorMap: Map[Identifier, Rational], + //trackInitialErrs: Boolean, + trackRoundoffErrs: Boolean, + uniformPrecision: Precision): (Interval, Rational) = { + + // TODO: this should not be hardcoded + val numSplits = 10 + + val inputsSubdiv: Seq[Map[Identifier, Interval]] = initValMap.foldLeft(Seq(Map[Identifier, Interval]()))({ + case (currSeq: Seq[Map[Identifier, Interval]], (id, intrvl)) => + val xlo = intrvl.xlo + val splitWidth = intrvl.width / numSplits + val splits: Seq[Interval] = (0 until numSplits).map(i => + Interval(xlo + i * splitWidth, xlo + (i+1) * splitWidth)) + + currSeq.flatMap( m => + splits.map(i => m + (id -> i)) + ) + }) + + val allIDs = inputParams + val missingIDs = allIDs -- initErrorMap.keySet + + // the initErrorMap will simply be empty, if initial errors are not being tracked + def getInputErrors(rangeMap: Map[Identifier, Interval]): Map[Identifier, Rational] = + if (trackRoundoffErrs){ + + initErrorMap ++ missingIDs.map(id => (id -> uniformPrecision.absRoundoff(rangeMap(id)))) + + } else { + + initErrorMap ++ missingIDs.map( id => (id -> Rational.zero)) + + } + + + val (resError, resInterval,_ ,_) = uniformRoundoff_IA_AA(expr, + inputsSubdiv.head, getInputErrors(inputsSubdiv.head), uniformPrecision) + // val (resInterval, resError) = evaluate[Interval, AffineForm]( + // expr, inputsSubdiv.head, + // getInputErrors(inputsSubdiv.head).map(x => (x._1 -> AffineForm.fromError(x._2))), + // Interval.apply, AffineForm.zero, AffineForm.fromError, AffineForm.apply, + // false, trackRoundoffErrs, uniformPrecision) + + var maxInterval = resInterval + var maxError = resError + + for (m <- inputsSubdiv.tail) { + + val (tmpError, tmpRange,_, _) = uniformRoundoff_IA_AA(expr, m, + getInputErrors(m), uniformPrecision) + + // evaluate[Interval, AffineForm]( + // expr, m, + // getInputErrors(m).map(x => (x._1 -> AffineForm.fromError(x._2))), + // Interval.apply, AffineForm.zero, AffineForm.fromError, AffineForm.apply, + // false, trackRoundoffErrs, uniformPrecision) + + maxInterval = maxInterval.union(tmpRange) + maxError = max(maxError, tmpError) + } + (maxInterval, maxError) + } + + + /** + * Uses interval subdivision to compute ranges and affine arithmetic for errors. + * + * Unlike the other subdivision method, this one annotates the trees with + * information. For this, ranges are computed first with interval subdivision, + * and errors are computed in a second step. Due to this separation, the computed + * errors may not be the tightest possible. + * This methods is suitable for fixed-point arithmetic. + * + * TODO: This needs to be updated to work without annotating trees. + * + * @param expr expression to be analysed + * @param inputRanges initial ranges of free variables + * @param varErrorMap input errors (incl. roundoff) of free variables + * @param trackRoundoffErrs whether or not to track roundoff errors or to + * propagate only input errors (as given in varErrorMap) + * @param uniformPrecision the default precision to use when doing absolute + * errors calculations + */ + // def evaluateSubdiv( + // expr: Expr, + // inputRanges: Map[Identifier, Interval], + // varErrorMap: Map[Identifier, Rational], + // trackRoundoffErrs: Boolean, + // uniformPrecision: Precision): (Interval, Rational) = { + + // @inline + // def computeNewError(range: Interval, propagatedError: AffineForm): AffineForm = + // if (trackRoundoffErrs) { + // val actualRange: Interval = range + propagatedError.toInterval + // val rndoff = uniformPrecision.absRoundoff(actualRange) + // propagatedError +/- rndoff + // } else { + // propagatedError + // } + + // val numSplits = 10 + + // //val inputRanges: Map[Identifier, Interval] = context.inputRanges(fncId) + + // val inputsSubdiv: Seq[Map[Identifier, Interval]] = inputRanges.foldLeft(Seq(Map[Identifier, Interval]()))({ + // case (currSeq: Seq[Map[Identifier, Interval]], (id, intrvl)) => + // val xlo = intrvl.xlo + // val splitWidth = intrvl.width / numSplits + // val splits: Seq[Interval] = (0 until numSplits).map(i => + // Interval(xlo + i * splitWidth, xlo + (i+1) * splitWidth)) + + // currSeq.flatMap( m => + // splits.map(i => m + (id -> i)) + // ) + // }) + + // // map of maximum seen ranges of all intermediate expressions + // var currentRanges: Map[Expr, Interval] = Map() + + // // there is probably a helper function for this + // def attachInterval(e: Expr): Unit = (e: @unchecked) match { + + // case x @ Variable(id) => + // if (!x.hasInterval) { + // x.interval = currentRanges(x) + // } + // return + + // case x @ RealLiteral(r) => + // x.interval = Interval(r) + // return + + // case x @ Plus(lhs, rhs) => + // attachInterval(lhs) + // attachInterval(rhs) + // x.interval = currentRanges(x) + + // case x @ Minus(lhs, rhs) => + // attachInterval(lhs) + // attachInterval(rhs) + // x.interval = currentRanges(x) + + // case x @ Times(lhs, rhs) => + // attachInterval(lhs) + // attachInterval(rhs) + // x.interval = currentRanges(x) + + // case x @ Division(lhs, rhs) => + // attachInterval(lhs) + // attachInterval(rhs) + // x.interval = currentRanges(x) + + // case x @ UMinus(t) => + // attachInterval(t) + // x.interval = currentRanges(x) + + // case x @ Sqrt(t) => + // attachInterval(t) + // x.interval = currentRanges(x) + + // // TODO: we may be missing intervals on the let-defined vars + // case x @ Let(id, value, body) => + // attachInterval(value) + // attachInterval(body) + // // the range of the let is the range of the final result returned + // x.interval = currentRanges(body) + // currentRanges = currentRanges + (x -> x.interval) + // } + + // // assumes that the ranges are attached + // def computeAndAttachError(e: Expr, errorMap: Map[Identifier, AffineForm]): AffineForm = (e: @unchecked) match { + // case x @ Variable(id) => + // val err = errorMap(id) + // x.absError = maxAbs(err.toInterval) + // err + + // case x @ RealLiteral(r) => + // val err = if (isExactInFloats(r, uniformPrecision)) { + // AffineForm.zero + // } else { + // AffineForm.fromError(uniformPrecision.absRoundoff(r)) + // } + // x.absError = maxAbs(err.toInterval) + // err + + // case x @ Plus(lhs, rhs) => + // val lErr = computeAndAttachError(lhs, errorMap) + // val rErr = computeAndAttachError(rhs, errorMap) + + // val propagatedError = lErr + rErr + // val newError = computeNewError(x.interval, propagatedError) + + // x.absError = maxAbs(newError.toInterval) + // newError + + // case x @ Minus(lhs, rhs) => + // val lErr = computeAndAttachError(lhs, errorMap) + // val rErr = computeAndAttachError(rhs, errorMap) + + // val propagatedError = lErr + rErr + // val newError = computeNewError(x.interval, propagatedError) + + // x.absError = maxAbs(newError.toInterval) + // newError + + // case x @ Times(lhs, rhs) => + // // propagated error + // val lAA = AffineForm(lhs.asInstanceOf[NumAnnotation].interval) + // val rAA = AffineForm(rhs.asInstanceOf[NumAnnotation].interval) + // val lErr = computeAndAttachError(lhs, errorMap) + // val rErr = computeAndAttachError(rhs, errorMap) + + // var propagatedError = lAA*rErr + rAA*lErr + lErr*rErr + + // // new error + // val newError = computeNewError(x.interval, propagatedError) + + // x.absError = maxAbs(newError.toInterval) + // newError + + // case x @ Division(lhs, rhs) => + // // propagated error + + // val rhsInterval = rhs.asInstanceOf[NumAnnotation].interval + // val rErr = computeAndAttachError(rhs, errorMap) //yErr + + // val rInt = rhsInterval + rErr.toInterval // the actual interval, incl errors + // val a = min(abs(rInt.xlo), abs(rInt.xhi)) + // val errorMultiplier = -one / (a*a) + + // val invErr = rErr * AffineForm(errorMultiplier) + // val lAA = AffineForm(lhs.asInstanceOf[NumAnnotation].interval) + // val inverse: Interval = rhsInterval.inverse + // val kAA = AffineForm(inverse) + + // val lErr = computeAndAttachError(lhs, errorMap) //xErr + + // // multiplication with inverse + // var propagatedError = lAA*invErr + kAA*lErr + lErr*invErr + + // // new error + // val newError = computeNewError(x.interval, propagatedError) + + // x.absError = maxAbs(newError.toInterval) + // newError + + // case x @ UMinus(t) => + // val newError = - computeAndAttachError(t, errorMap) // no additional roundoff error + // x.absError = maxAbs(newError.toInterval) + // newError + + // case x @ Sqrt(t) => + + // // TODO: Not supported for fixed-points, add exception + // //case FPPrecision(_) => throw UnsupportedRealFragmentException("Sqrt not supported for fixed-points.") + + // val tInterval = t.asInstanceOf[NumAnnotation].interval + // val tError = computeAndAttachError(t, errorMap) + + // // propagated existing errors + // val a = min(abs(tInterval.xlo), abs(tInterval.xhi)) + // val errorMultiplier = Rational(1L, 2L) / sqrtDown(a) + // val propagatedError = tError * AffineForm(errorMultiplier) + + // // new roundoff + // val newError = computeNewError(x.interval, propagatedError) + + // x.absError = maxAbs(newError.toInterval) + // newError + + + // case Let(id, value, body) => + // val aform = computeAndAttachError(value, errorMap) + // val newErrorMap = errorMap + (id -> aform) + // computeAndAttachError(body, newErrorMap) + // } + + // // collect possible intervals + // // first intervals need to populate the map + // evalInterval(expr, inputsSubdiv.head, + // (x: Expr, i: Interval) => currentRanges = currentRanges + (x -> i)) + + // // the next intervals update the map + // inputsSubdiv.tail.foreach( subdiv => + // evalInterval(expr, subdiv, + // (x: Expr, i: Interval) => currentRanges = currentRanges + (x -> currentRanges(x).union(i)))) + + // currentRanges = currentRanges ++ inputRanges.map(x => (Variable(x._1) -> x._2)) + + // // traverse the tree once more and attach the correct ranges + // attachInterval(expr) + + // // Now that we have the ranges, we can compute the errors + // val resError = computeAndAttachError(expr, varErrorMap.map(x => (x._1 -> AffineForm.fromError(x._2)))) + + // (currentRanges(expr), maxAbs(resError.toInterval)) + // } + +} diff --git a/src/main/scala/daisy/utils/RangeEvaluators.scala b/src/main/scala/daisy/utils/RangeEvaluators.scala new file mode 100644 index 0000000000000000000000000000000000000000..4fa1516221868d780c54c6b0727d14aee8510625 --- /dev/null +++ b/src/main/scala/daisy/utils/RangeEvaluators.scala @@ -0,0 +1,94 @@ + +package daisy +package utils + +import scala.collection.immutable.Seq + + +import lang.Trees._ +import lang.Identifiers._ + +trait RangeEvaluators { + + // TODO: this function probably subsumes other evaluators, probably remove them? + // assert stmt that precision map cannot mix floats and fps + // regression tests for mixed-precision and rewriting (separately) + // TODO: fix fixed-point arithmetic code generation + + /** + * Evaluates the range of this expression using the given RangeArithmetic + * and saves the intermediate ranges in a map. + * + * TODO: documentation + * + * @param {[type]} Rational [description] + * @return {[type]} [description] + */ + def evalRange[T <: RangeArithmetic[T]]( + expr: Expr, + initValMap: Map[Identifier, T], + rangeFromReal: Rational => T): (T, Map[Expr, T]) = { + + // TODO: we could even do easy caching now, since we can just check whether the + // stuff is already computed... + var intermediateRanges: Map[Expr, T] = Map.empty + + def eval(e: Expr, valMap: Map[Identifier, T]): T = (e: @unchecked) match { + + case x @ RealLiteral(r) => + val range = rangeFromReal(r) + intermediateRanges += (x -> range) + range + + case x @ Variable(id) => + // we update the map each time we encounter a variable, in principle + // we only need to do this once though. Perhaps optimise? + val range = valMap(id) + intermediateRanges += (x -> range) + range + + case x @ Plus(lhs, rhs) => + val range = eval(lhs, valMap) + eval(rhs, valMap) + intermediateRanges += (x -> range) + range + + case x @ Minus(lhs, rhs) => + val range = eval(lhs, valMap) - eval(rhs, valMap) + intermediateRanges += (x -> range) + range + + case x @ Times(lhs, rhs) => + val range = eval(lhs, valMap) * eval(rhs, valMap) + intermediateRanges += (x -> range) + range + + + case x @ Division(lhs, rhs) => + val range = eval(lhs, valMap) / eval(rhs, valMap) + intermediateRanges += (x -> range) + range + + case x @ UMinus(t) => + val range = - eval(t, valMap) + intermediateRanges += (x -> range) + range + + case x @ Sqrt(t) => + val range = eval(t, valMap).squareRoot + intermediateRanges += (x -> range) + range + + + case x @ Let(id, value, body) => + val valueRange = eval(value, valMap) + // TODO: do we need the ranges also of the Let's? + eval(body, valMap + (id -> valueRange)) + + } + val res = eval(expr, initValMap) + (res, intermediateRanges) + } + + + +} \ No newline at end of file diff --git a/src/main/scala/daisy/utils/Rational.scala b/src/main/scala/daisy/utils/Rational.scala index 19648e57e82d45c7c03925272ad3076e8c0e5c56..0bcd42c75e4f570b61da908b52bd2f6d689d6b40 100755 --- a/src/main/scala/daisy/utils/Rational.scala +++ b/src/main/scala/daisy/utils/Rational.scala @@ -485,6 +485,10 @@ class Rational private(val n: BigInt, val d: BigInt) extends ScalaNumber with Sc case _ => false } + override def hashCode(): Int = { + this.doubleValue.hashCode + } + override def byteValue(): Byte = Predef.double2Double(doubleValue).byteValue override def doubleValue(): Double = { val bigN = new java.math.BigDecimal(n.bigInteger, mathContext) diff --git a/src/main/scala/daisy/utils/RoundoffEvaluators.scala b/src/main/scala/daisy/utils/RoundoffEvaluators.scala new file mode 100644 index 0000000000000000000000000000000000000000..5b58909541520d43a2689ce2ebbe1f59101c5b55 --- /dev/null +++ b/src/main/scala/daisy/utils/RoundoffEvaluators.scala @@ -0,0 +1,316 @@ + +package daisy +package utils + +import scala.collection.immutable.Seq + +import lang.Trees._ +import lang.Identifiers._ +import lang.TreeOps.allVariablesOf +import FinitePrecision._ +import Rational.{min, abs, sqrtDown, one} + +trait RoundoffEvaluators extends RangeEvaluators { + + /** + * Calculate static error using interval ranges and interval errors + * + * @param expr Expression whose erros is to be calculated + * @param inputValMap Map from function identifier to its input value interval + * @param inputErrorMap Map from function + */ + def uniformRoundoff_IA_IA( + expr: Expr, + inputValMap: Map[Identifier, Interval], + inputErrorMap: Map[Identifier, Rational], + uniformPrecision: Precision, + trackRoundoffErrors: Boolean = true): (Rational, Interval, Map[Expr, Interval], Map[Expr, Interval]) = { + + val (resRange, intermediateRanges) = evalRange[Interval](expr, + inputValMap.map(x => (x._1 -> Interval(x._2))), Interval.apply) + + val (resRoundoff, allErrors) = evalRoundoff[Interval](expr, intermediateRanges, + allVariablesOf(expr).map(id => (id -> uniformPrecision)).toMap, + inputErrorMap.map(x => (x._1 -> Interval.fromError(x._2))), + zeroError = Interval.zero, + fromError = Interval.fromError, + interval2T = Interval.apply, + constantsPrecision = uniformPrecision, + trackRoundoffErrors) + + (Interval.maxAbs(resRoundoff.toInterval),resRange, intermediateRanges, allErrors) + } + + /** + * Calculates the roundoff error for a given uniform precision + * using interval arithmetic for ranges and affine arithmetic for errors. + * + * @param expr expression for which to compute roundoff + * @param inputValMap real-valued ranges of all input variables + * @param inputErrorMap errors of all input variables (incl. roundoff) + * @param uniformPrecision precision for the entire computation + * + * @return (max. absolute roundoff error bound, real-valued result interval) + */ + def uniformRoundoff_IA_AA( + expr: Expr, + inputValMap: Map[Identifier, Interval], + inputErrorMap: Map[Identifier, Rational], + uniformPrecision: Precision, + trackRoundoffErrors: Boolean = true): (Rational, Interval, Map[Expr, Interval], Map[Expr, Interval]) = { + + val (resRange, intermediateRanges) = evalRange[Interval](expr, inputValMap, Interval.apply) + //println(intermediateRanges.mkString("\n")) + + val (resRoundoff, allErrors) = evalRoundoff[AffineForm](expr, intermediateRanges, + allVariablesOf(expr).map(id => (id -> uniformPrecision)).toMap, + inputErrorMap.map(x => (x._1 -> AffineForm.fromError(x._2))), + zeroError = AffineForm.zero, + fromError = AffineForm.fromError, + interval2T = AffineForm.apply, + constantsPrecision = uniformPrecision, + trackRoundoffErrors) + + (Interval.maxAbs(resRoundoff.toInterval), resRange, intermediateRanges, + allErrors.map(x => (x._1 -> x._2.toInterval))) + } + + /** + * Calculates the roundoff error for a given uniform precision + * using affine arithmetic for ranges and affine arithmetic for errors. + * + * @param expr expression for which to compute roundoff + * @param inputValMap real-valued ranges of all input variables + * @param inputErrorMap errors of all input variables (incl. roundoff) + * @param uniformPrecision precision for the entire computation + */ + def uniformRoundoff_AA_AA( + expr: Expr, + inputValMap: Map[Identifier, Interval], + inputErrorMap: Map[Identifier, Rational], + uniformPrecision: Precision, + trackRoundoffErrors: Boolean = true): (Rational, Interval, Map[Expr, Interval], Map[Expr, Interval]) = { + + val (resRange, intermediateRanges) = evalRange[AffineForm](expr, + inputValMap.map(x => (x._1 -> AffineForm(x._2))), AffineForm.apply) + + val (resRoundoff, allErrors) = evalRoundoff[AffineForm](expr, + intermediateRanges.map(x => (x._1 -> x._2.toInterval)), + allVariablesOf(expr).map(id => (id -> uniformPrecision)).toMap, + inputErrorMap.map(x => (x._1 -> AffineForm.fromError(x._2))), + zeroError = AffineForm.zero, + fromError = AffineForm.fromError, + interval2T = AffineForm.apply, + constantsPrecision = uniformPrecision, + trackRoundoffErrors) + + (Interval.maxAbs(resRoundoff.toInterval), resRange.toInterval, + intermediateRanges.map(x => (x._1 -> x._2.toInterval)), + allErrors.map(x => (x._1 -> x._2.toInterval))) + } + + /** + * Calculates the roundoff error for a given uniform precision + * using SMTRange for ranges and affine arithmetic for errors. + * + * @param expr expression for which to compute roundoff + * @param inputValMap real-valued ranges of all input variables + * @param inputErrorMap errors of all input variables (incl. roundoff) + * @param uniformPrecision precision for the entire computation + */ + def uniformRoundoff_SMT_AA( + expr: Expr, + inputValMap: Map[Identifier, Interval], + inputErrorMap: Map[Identifier, Rational], + uniformPrecision: Precision, + trackRoundoffErrors: Boolean = true): (Rational, Interval, Map[Expr, Interval], Map[Expr, Interval]) = { + + val (resRange, intermediateRanges) = evalRange[SMTRange](expr, + inputValMap.map({ case (id, int) => (id -> SMTRange(Variable(id), int)) }), + SMTRange.apply) + + val (resRoundoff, allErrors) = evalRoundoff[AffineForm](expr, + intermediateRanges.map(x => (x._1 -> x._2.toInterval)), + allVariablesOf(expr).map(id => (id -> uniformPrecision)).toMap, + inputErrorMap.map(x => (x._1 -> AffineForm.fromError(x._2))), + zeroError = AffineForm.zero, + fromError = AffineForm.fromError, + interval2T = AffineForm.apply, + constantsPrecision = uniformPrecision, + trackRoundoffErrors) + + (Interval.maxAbs(resRoundoff.toInterval), resRange.toInterval, + intermediateRanges.map(x => (x._1 -> x._2.toInterval)), + allErrors.map(x => (x._1 -> x._2.toInterval))) + } + + /** + Computes the absolute roundoff error for the given expression. + + The ranges of all the intermediate expressions have to be given in rangeMap. + Allows mixed-precision by providing (possibly different) precisions for + all declared variables (input parameters as well as locally defined variables.) + Constants are assumed to be all in one precision, given by the user. + + */ + def evalRoundoff[T <: RangeArithmetic[T]]( + expr: Expr, + rangeMap: Map[Expr, Interval], + precisionMap: Map[Identifier, Precision], + initErrorMap: Map[Identifier, T], + zeroError: T, + fromError: Rational => T, + interval2T: Interval => T, + constantsPrecision: Precision, + trackRoundoffErrors: Boolean // if false, propagate only initial errors + ): (T, Map[Expr, T]) = { + + + var intermediateErrors: Map[Expr, T] = Map.empty + + // TODO: check the effectiveness of this + //@inline + def computeNewError(range: Interval, propagatedError: T, prec: Precision): T = + if (trackRoundoffErrors) { + val actualRange: Interval = range + propagatedError.toInterval + val rndoff = prec.absRoundoff(actualRange) + propagatedError +/- rndoff + } else { + propagatedError + } + + def eval(e: Expr, errorMap: Map[Identifier, T]): (T, Precision) = (e: @unchecked) match { + + case x @ RealLiteral(r) => + //val rndoff = if (isExactInFloats(r, constantsPrecision) || !trackRoundoffErrors) { + val rndoff = if (!trackRoundoffErrors) { + zeroError + } else { + fromError(constantsPrecision.absRoundoff(r)) + } + intermediateErrors += (x -> rndoff) + (rndoff, constantsPrecision) + + case x @ Variable(id) => + // TODO: if the error is just a roundoff, then we can also compute it here... + val rndoff = errorMap(id) + intermediateErrors += (x -> rndoff) + (rndoff, precisionMap(id)) + + case x @ Plus(lhs, rhs) => + val range = rangeMap(x) + + val (rndoffLhs, precLhs) = eval(lhs, errorMap) + val (rndoffRhs, precRhs) = eval(rhs, errorMap) + val propagatedError = rndoffLhs + rndoffRhs + + val prec = getUpperBound(precLhs, precRhs) // Scala semantics + val rndoff = computeNewError(range, propagatedError, prec) + intermediateErrors += (x -> rndoff) + (rndoff, prec) + + case x @ Minus(lhs, rhs) => + val range = rangeMap(x) + + val (rndoffLhs, precLhs) = eval(lhs, errorMap) + val (rndoffRhs, precRhs) = eval(rhs, errorMap) + val propagatedError = rndoffLhs - rndoffRhs + + val prec = getUpperBound(precLhs, precRhs) + val rndoff = computeNewError(range, propagatedError, prec) + intermediateErrors += (x -> rndoff) + (rndoff, prec) + + case x @ Times(lhs, rhs) => + val range = rangeMap(x) + val rangeLhs = rangeMap(lhs) + val rangeRhs = rangeMap(rhs) + + val (rndoffLhs, precLhs) = eval(lhs, errorMap) + val (rndoffRhs, precRhs) = eval(rhs, errorMap) + + val propagatedError = interval2T(rangeLhs) * rndoffRhs + + interval2T(rangeRhs) * rndoffLhs + + rndoffLhs * rndoffRhs + + val prec = getUpperBound(precLhs, precRhs) + val rndoff = computeNewError(range, propagatedError, prec) + intermediateErrors += (x -> rndoff) + (rndoff, prec) + + + case x @ Division(lhs, rhs) => + val range = rangeMap(x) + val rangeLhs = rangeMap(lhs) + val rangeRhs = rangeMap(rhs) + + val (rndoffLhs, precLhs) = eval(lhs, errorMap) + val (rndoffRhs, precRhs) = eval(rhs, errorMap) + + // inverse, i.e. we are computing x * (1/y) + val rightInterval = rangeRhs + rndoffRhs.toInterval // the actual interval, incl errors + + //the actual error interval can now contain 0, check this + if (rightInterval.xlo <= 0 && rightInterval.xhi >= 0) + throw DivisionByZeroException("trying to divide by error interval containing 0") + + val a = min(abs(rightInterval.xlo), abs(rightInterval.xhi)) + val errorMultiplier: Rational = -one / (a*a) + val invErr = rndoffRhs * errorMultiplier + + // error propagation + val inverse: Interval = rangeRhs.inverse + + var propagatedError = interval2T(rangeLhs) * invErr + + interval2T(inverse) * rndoffLhs + + rndoffLhs * invErr + + val prec = getUpperBound(precLhs, precRhs) + val rndoff = computeNewError(range, propagatedError, prec) + intermediateErrors += (x -> rndoff) + (rndoff, prec) + + case x @ UMinus(t) => + val (rndoff, prec) = eval(t, errorMap) + intermediateErrors += (x -> - rndoff) + (- rndoff, prec) + + case x @ Sqrt(t) => + // TODO: needs to fail for fixed-point precision + val range = rangeMap(x) + val rangeT = rangeMap(t) + val (errorT, prec) = eval(t, errorMap) + + val tInterval = rangeT + val a = min(abs(tInterval.xlo), abs(tInterval.xhi)) + val errorMultiplier = Rational(1L, 2L) / sqrtDown(a) + val propagatedError = errorT * errorMultiplier + + // TODO: check that this operation exists for this precision + val rndoff = computeNewError(range, propagatedError, prec) + + intermediateErrors += (x -> rndoff) + (rndoff, prec) + + + case x @ Let(id, value, body) => + val (valueRndoff, valuePrec) = eval(value, errorMap) + + // downcast required + val idPrec = precisionMap(id) + val rndoff = if (idPrec < valuePrec) { // we need to cast down + val valueRange = rangeMap(value) + computeNewError(valueRange, valueRndoff, idPrec) + } else { + valueRndoff + } + + eval(body, errorMap + (id -> rndoff)) + + } + val (resRndoff, resPrecision) = eval(expr, initErrorMap) + (resRndoff, intermediateErrors) + } + + +} diff --git a/src/test/resources/AbsErrorRegressionFunctions.scala b/src/test/resources/AbsErrorRegressionFunctions.scala index e8577af97cd0a209dd230056267b1bbfd55d9897..155c7c4a48ccb93f4a3777a647010094b955e04a 100644 --- a/src/test/resources/AbsErrorRegressionFunctions.scala +++ b/src/test/resources/AbsErrorRegressionFunctions.scala @@ -106,6 +106,41 @@ object RangeRegressionFunctions { } + def kepler0(x1: Real, x2: Real, x3: Real, x4: Real, x5: Real, x6: Real): Real = { + require(4 <= x1 && x1 <= 6.36 && 4 <= x2 && x2 <= 6.36 && 4 <= x3 && x3 <= 6.36 && + 4 <= x4 && x4 <= 6.36 && 4 <= x5 && x5 <= 6.36 && 4 <= x6 && x6 <= 6.36) + + x2 * x5 + x3 * x6 - x2 * x3 - x5 * x6 + x1 * (-x1 + x2 + x3 - x4 + x5 + x6) + + } // 1.15e-15 + + + def kepler1(x1: Real, x2: Real, x3: Real, x4: Real): Real = { + require(4 <= x1 && x1 <= 6.36 && 4 <= x2 && x2 <= 6.36 && 4 <= x3 && x3 <= 6.36 && + 4 <= x4 && x4 <= 6.36) + + x1 * x4 * (-x1 + x2 + x3 - x4) + x2 * (x1 - x2 + x3 + x4) + x3 * (x1 + x2 - x3 + x4) - + x2 * x3 * x4 - x1 * x3 - x1 * x2 - x4 + + } // 4.50e–13 + + def kepler2(x1: Real, x2: Real, x3: Real, x4: Real, x5: Real, x6: Real): Real = { + require(4 <= x1 && x1 <= 6.36 && 4 <= x2 && x2 <= 6.36 && 4 <= x3 && x3 <= 6.36 && + 4 <= x4 && x4 <= 6.36 && 4 <= x5 && x5 <= 6.36 && 4 <= x6 && x6 <= 6.36) + + x1 * x4 * (-x1 + x2 + x3 - x4 + x5 + x6) + x2 * x5 * (x1 - x2 + x3 + x4 - x5 + x6) + + x3* x6 * (x1 + x2 - x3 + x4 + x5 - x6) - x2 * x3 * x4 - + x1* x3* x5 - x1 * x2 * x6 - x4 * x5 * x6 + + } //2.08e–12 + + def himmilbeau(x1: Real, x2: Real) = { + require(-5 <= x1 && x1 <= 5 && -5 <= x2 && x2 <= 5) + + (x1*x1 + x2 - 11)*(x1 * x1 + x2 - 11) + (x1 + x2*x2 - 7)*(x1 + x2*x2 - 7) + + } //1.43e–12 + /*def jetEngine(x1: Real, x2: Real) = { require(-5 <= x1 && x1 <= 5 && -20 <= x2 && x2 <= 5) x1 + ( diff --git a/src/test/scala/regression/AbsErrorAnalysisRegressionTest.scala b/src/test/scala/regression/AbsErrorAnalysisRegressionTest.scala index 2b86a943cef3a0b833dd9798e8c1a9e1e64bf05c..f8dd43fdfdf0b0fb43f4bc67dddcaecc02957a1f 100644 --- a/src/test/scala/regression/AbsErrorAnalysisRegressionTest.scala +++ b/src/test/scala/regression/AbsErrorAnalysisRegressionTest.scala @@ -7,6 +7,23 @@ import org.scalatest.FunSuite import lang.Trees.{Expr, Let} import lang.NumAnnotation +/** + Regression test for the basic absolute error and range computations. + + Times taken on Eva's Linux desktop (before refactoring | after refactoring): + int-affine-fixed32 425 413 438 | 361 365 + int-affine-float32 307 304 321 | 281 278 + int-affine-float64 593 590 603 | 505 509 + + affine-affine-fixed32 453 434 461 | 382 382 + affine-affine-float32 329 334 347 | 288 285 + affine-affine-float64 625 635 650 | 527 522 + + smt-affine-fixed32 47078 46967 47108 | 50006 49934 + smt-affine-float32 46235 46350 46273 | 49413 49287 + smt-affine-float64 46622 46559 46486 | 49481 49474 + +*/ class AbsErrorAnalysisRegressionTest extends FunSuite { val inputFile: String = "src/test/resources/AbsErrorRegressionFunctions.scala" @@ -17,13 +34,58 @@ class AbsErrorAnalysisRegressionTest extends FunSuite { files = List(inputFile) ) - // default analysis: IA - AA val inputPrg = frontend.ExtractionPhase(context) val pipeline = analysis.SpecsProcessingPhase andThen analysis.RangeErrorPhase // Regression results // fnc.id -> (abs error, range) - val intAA: Map[String, (String, String)] = Map( + + val intAA_Fixed32: Map[String, (String, String)] = Map( + ("doppler", ("1.75791622644632e-06", "[-158.7191444098274, -0.02944244059231351]")), + ("sine", ("5.172930285837377e-09", "[-2.3011348046703466, 2.3011348046703466]")), + ("sineOrder3", ("7.483225706087375e-09", "[-2.9419084189651277, 2.9419084189651277]")), + ("sqroot", ("4.015397284953641e-06", "[-402.125, 68.5]")), + ("bspline0", ("8.537123605701497e-10", "[0.0, 0.16666666666666666]")), + ("bspline1", ("3.337239226446336e-09", "[-0.3333333333333333, 1.1666666666666667]")), + ("bspline2", ("3.104408582665976e-09", "[-0.3333333333333333, 1.1666666666666667]")), + ("bspline3", ("6.208817165548793e-10", "[-0.16666666666666666, 0.0]")), + ("rigidBody1", ("1.3485550881126018e-06", "[-705.0, 705.0]")), + ("rigidBody2", ("0.0001531802118132919", "[-58740.0, 58740.0]")), + ("verhulst", ("3.009066920607247e-09", "[0.3148936170212766, 1.1008264462809918]")), + ("predatorPrey", ("1.647071564900709e-09", "[0.03727705922396188, 0.35710168263424846]")), + ("carbonGas", ("17.69469347003061", "[2097409.2, 3.434323E7]")), + ("turbine1", ("5.059999504388075e-07", "[-58.32912689020381, -1.5505285721480735]")), + ("turbine2", ("6.241048750415222e-07", "[-29.43698909090909, 80.993]")), + ("turbine3", ("3.6766172442271937e-07", "[0.4660958448753463, 40.375126890203816]")), + ("kepler0", ("4.379078746518522e-07", "[-35.7792, 159.8176]")), + ("kepler1", ("2.021807432895287e-06", "[-490.320768, 282.739712]")), + ("kepler2", ("1.037595868540202e-05", "[-871.597824, 1860.323072]")), + ("himmilbeau", ("9.75281000696604e-06", "[-1630.0, 3050.0]")) + ) + + val intAA_Float32: Map[String, (String, String)] = Map( + ("doppler", ("0.0002250138141196529", "[-158.7191444098274, -0.02944244059231351]")), + ("sine", ("6.064885588214318e-07", "[-2.3011348046703466, 2.3011348046703466]")), + ("sineOrder3", ("7.790389913381473e-07", "[-2.9419084189651277, 2.9419084189651277]")), + ("sqroot", ("0.00016796589915202306", "[-402.125, 68.5]")), + ("bspline0", ("8.69234485871102e-08", "[0.0, 0.16666666666666666]")), + ("bspline1", ("4.2716663302873797e-07", "[-0.3333333333333333, 1.1666666666666667]")), + ("bspline2", ("3.973643085686263e-07", "[-0.3333333333333333, 1.1666666666666667]")), + ("bspline3", ("5.712112027822517e-08", "[-0.16666666666666666, 0.0]")), + ("rigidBody1", ("0.00017261505240639963", "[-705.0, 705.0]")), + ("rigidBody2", ("0.01960706761337861", "[-58740.0, 58740.0]")), + ("verhulst", ("3.199529347654797e-07", "[0.3148936170212766, 1.1008264462809918]")), + ("predatorPrey", ("1.3291118406629398e-07", "[0.03727705922396188, 0.35710168263424846]")), + ("carbonGas", ("22.995233482473232", "[2097409.2, 3.434323E7]")), + ("turbine1", ("5.095029323813888e-05", "[-58.32912689020381, -1.5505285721480735]")), + ("turbine2", ("7.46823932469574e-05", "[-29.43698909090909, 80.993]")), + ("turbine3", ("3.797059491586537e-05", "[0.4660958448753463, 40.375126890203816]")), + ("kepler0", ("5.605220905522401e-05", "[-35.7792, 159.8176]")), + ("kepler1", ("0.0002587913631247376", "[-490.320768, 282.739712]")), + ("kepler2", ("0.0013281227815583414", "[-871.597824, 1860.323072]")), + ("himmilbeau", ("0.0012483597718073725", "[-1630.0, 3050.0]")) + ) + val intAA_Float64: Map[String, (String, String)] = Map( ("doppler", ("4.1911988101104756e-13", "[-158.7191444098274, -0.02944244059231351]")), ("sine", ("1.1296729607621835e-15", "[-2.3011348046703466, 2.3011348046703466]")), ("sineOrder3", ("1.4510731312549944e-15", "[-2.9419084189651277, 2.9419084189651277]")), @@ -39,10 +101,61 @@ class AbsErrorAnalysisRegressionTest extends FunSuite { ("carbonGas", ("4.2831938128927477e-08", "[2097409.2, 3.434323E7]")), ("turbine1", ("9.490226544267962e-14", "[-58.32912689020381, -1.5505285721480735]")), ("turbine2", ("1.3910672706969435e-13", "[-29.43698909090909, 80.993]")), - ("turbine3", ("7.072570725135768e-14", "[0.4660958448753463, 40.375126890203816]")) + ("turbine3", ("7.072570725135768e-14", "[0.4660958448753463, 40.375126890203816]")), + ("kepler0", ("1.0440537323574973e-13", "[-35.7792, 159.8176]")), + ("kepler1", ("4.820364551960666e-13", "[-490.320768, 282.739712]")), + ("kepler2", ("2.4738213255659507e-12", "[-871.597824, 1860.323072]")), + ("himmilbeau", ("2.3252511027749283e-12", "[-1630.0, 3050.0]")) + ) - val aaAA: Map[String, (String, String)] = Map( + val aaAA_Fixed32: Map[String, (String, String)] = Map( + ("doppler", ("5.390735253968273e-06", "[-282.4762120545625, 133.53913086064713]")), + ("sine", ("3.5841039195923494e-09", "[-1.6540002686363373, 1.6540002686363373]")), + ("sineOrder3", ("6.536979105441602e-09", "[-1.909859317102744, 1.909859317102744]")), + ("sqroot", ("4.004454244701909e-06", "[-365.875, 195.04296875]")), + ("bspline0", ("8.537123605701497e-10", "[-0.0625, 0.16666666666666666]")), + ("bspline1", ("3.162616243705934e-09", "[-0.08333333333333333, 0.9791666666666666]")), + ("bspline2", ("2.8521753853742843e-09", "[-0.020833333333333332, 0.9166666666666666]")), + ("bspline3", ("6.208817165548793e-10", "[-0.16666666666666666, 0.0625]")), + ("rigidBody1", ("1.3485550881126018e-06", "[-705.0, 705.0]")), + ("rigidBody2", ("0.0001531802118132919", "[-57390.0, 58740.0]")), + ("verhulst", ("2.7762362769533775e-09", "[0.35217660784145, 0.9891620430819574]")), + ("predatorPrey", ("1.6642254365978505e-09", "[-0.0021029110350001693, 0.34533842311741914]")), + ("carbonGas", ("17.682124422310654", "[-6856426.768, 3.159230584E7]")), + ("turbine1", ("5.106469765682726e-07", "[-55.55054339020381, 38.36885388122605]")), + ("turbine2", ("6.326289514425209e-07", "[-80.18793076923077, 63.665174692307694]")), + ("turbine3", ("3.7230875055215864e-07", "[-27.227529265841437, 38.719203390203816]")), + ("kepler0", ("3.8016587502569254e-07", "[3.0664, 102.8708]")), + ("kepler1", ("1.842993498569115e-06", "[-288.15184, 89.927712]")), + ("kepler2", ("9.235736732009992e-06", "[-337.61856, 850.310096]")), + ("himmilbeau", ("4.194676880030224e-06", "[-775.0, 890.0]")) + ) + + val aaAA_Float32: Map[String, (String, String)] = Map( + ("doppler", ("0.0006900171033357679", "[-282.4762120545625, 133.53913086064713]")), + ("sine", ("4.031187839420683e-07", "[-1.6540002686363373, 1.6540002686363373]")), + ("sineOrder3", ("6.579194264554884e-07", "[-1.909859317102744, 1.909859317102744]")), + ("sqroot", ("0.00016656518999980138", "[-365.875, 195.04296875]")), + ("bspline0", ("8.69234485871102e-08", "[-0.0625, 0.16666666666666666]")), + ("bspline1", ("4.048148912379665e-07", "[-0.08333333333333333, 0.9791666666666666]")), + ("bspline2", ("3.6507845931528976e-07", "[-0.020833333333333332, 0.9166666666666666]")), + ("bspline3", ("5.712112027822517e-08", "[-0.16666666666666666, 0.0625]")), + ("rigidBody1", ("0.00017261505240639963", "[-705.0, 705.0]")), + ("rigidBody2", ("0.01960706761337861", "[-57390.0, 58740.0]")), + ("verhulst", ("2.901506123777844e-07", "[0.35217660784145, 0.9891620430819574]")), + ("predatorPrey", ("1.342975706563773e-07", "[-0.0021029110350001693, 0.34533842311741914]")), + ("carbonGas", ("21.386395374318933", "[-6856426.768, 3.159230584E7]")), + ("turbine1", ("5.1545112652219064e-05", "[-55.55054339020381, 38.36885388122605]")), + ("turbine2", ("7.577347502729098e-05", "[-80.18793076923077, 63.665174692307694]")), + ("turbine3", ("3.8565414329523206e-05", "[-27.227529265841437, 38.719203390203816]")), + ("kepler0", ("4.866123310307557e-05", "[3.0664, 102.8708]")), + ("kepler1", ("0.00023590317953098758", "[-288.15184, 89.927712]")), + ("kepler2", ("0.0011821743715241617", "[-337.61856, 850.310096]")), + ("himmilbeau", ("0.0005369187050519044", "[-775.0, 890.0]")) + ) + + val aaAA_Float64: Map[String, (String, String)] = Map( ("doppler", ("1.2852513956990104e-12", "[-282.4762120545625, 133.53913086064713]")), ("sine", ("7.508672360829453e-16", "[-1.6540002686363373, 1.6540002686363373]")), ("sineOrder3", ("1.2254703612493457e-15", "[-1.909859317102744, 1.909859317102744]")), @@ -58,87 +171,237 @@ class AbsErrorAnalysisRegressionTest extends FunSuite { ("carbonGas", ("3.983524363087597e-08", "[-6856426.768, 3.159230584E7]")), ("turbine1", ("9.601020280849505e-14", "[-55.55054339020381, 38.36885388122605]")), ("turbine2", ("1.4113902525335888e-13", "[-80.18793076923077, 63.665174692307694]")), - ("turbine3", ("7.183364461717313e-14", "[-27.227529265841437, 38.719203390203816]")) - ) + ("turbine3", ("7.183364461717313e-14", "[-27.227529265841437, 38.719203390203816]")), + ("kepler0", ("9.063860773039779e-14", "[3.0664, 102.8708]")), + ("kepler1", ("4.394038910504606e-13", "[-288.15184, 89.927712]")), + ("kepler2", ("2.201971227577815e-12", "[-337.61856, 850.310096]")), + ("himmilbeau", ("1.0000889005823412e-12", "[-775.0, 890.0]")) + ) + + + val smtAA_Fixed32: Map[String, (String, String)] = Map( + ("doppler", ("1.75791622644632e-06", "[-137.64316836703838, -0.02944244059231351]")), + ("sine", ("3.350232349298144e-09", "[-1.0043283753002443, 1.0043283753002443]")), + ("sineOrder3", ("6.536979105441602e-09", "[-1.0055351041384715, 1.0055351041384715]")), + ("sqroot", ("3.998633478610562e-06", "[-334.7934013614431, 1.6171710207127035]")), + ("bspline0", ("8.537123605701497e-10", "[0.0, 0.16666666666666666]")), + ("bspline1", ("2.949188153689887e-09", "[0.16617838541666666, 0.6671549479166666]")), + ("bspline2", ("2.599942188082592e-09", "[0.16666666666666666, 0.6668650309244791]")), + ("bspline3", ("6.208817165548793e-10", "[-0.16666666666666666, 0.0]")), + ("rigidBody1", ("1.3485550881126018e-06", "[-705.0, 705.0]")), + ("rigidBody2", ("0.0001531802118132919", "[-56049.1667175293, 58740.0]")), + ("verhulst", ("2.7762362769533775e-09", "[0.3640144188500088, 0.9473239405662036]")), + ("predatorPrey", ("1.647071564900709e-09", "[0.03727705922396188, 0.33711264367110555]")), + ("carbonGas", ("17.67125597003061", "[4301713.35625, 1.6740286809375E7]")), + ("turbine1", ("4.799229183495741e-07", "[-18.536291471169722, -1.986569412987437]")), + ("turbine2", ("5.570496496692077e-07", "[-28.564652726384942, 3.8296552477598733]")), + ("turbine3", ("3.3413411173656214e-07", "[0.5626821330545027, 11.434005458108444]")), + ("kepler0", ("3.8016587502569254e-07", "[20.802505677616598, 95.98217598131895]")), + ("kepler1", ("1.7833888537937244e-06", "[-278.40975267291554, -31.98171302891623]")), + ("kepler2", ("9.42228436899577e-06", "[-645.8671575, 1288.34451875]")), + ("himmilbeau", ("4.194676880030224e-06", "[-0.1416015625, 890.0]")) + ) - val smtAA: Map[String, (String, String)] = Map( + + val smtAA_Float32: Map[String, (String, String)] = Map( + ("doppler", ("0.0002250138141196529", "[-137.64316836703838, -0.02944244059231351]")), + ("sine", ("3.7318322294440997e-07", "[-1.0043283753002443, 1.0043283753002443]")), + ("sineOrder3", ("6.579194264554884e-07", "[-1.0055351041384715, 1.0055351041384715]")), + ("sqroot", ("0.000165820131940109", "[-334.7934013614431, 1.6171710207127035]")), + ("bspline0", ("8.69234485871102e-08", "[0.0, 0.16666666666666666]")), + ("bspline1", ("3.774960957159124e-07", "[0.16617838541666666, 0.6671549479166666]")), + ("bspline2", ("3.3030908319631186e-07", "[0.16666666666666666, 0.6668650309244791]")), + ("bspline3", ("5.712112027822517e-08", "[-0.16666666666666666, 0.0]")), + ("rigidBody1", ("0.00017261505240639963", "[-705.0, 705.0]")), + ("rigidBody2", ("0.01960706761337861", "[-56049.1667175293, 58740.0]")), + ("verhulst", ("2.901506123777844e-07", "[0.3640144188500088, 0.9473239405662036]")), + ("predatorPrey", ("1.3291118406629398e-07", "[0.03727705922396188, 0.33711264367110555]")), + ("carbonGas", ("19.995233482473232", "[4301713.35625, 1.6740286809375E7]")), + ("turbine1", ("4.7612433130717006e-05", "[-18.536291471169722, -1.986569412987437]")), + ("turbine2", ("6.609932439930116e-05", "[-28.564652726384942, 3.8296552477598733]")), + ("turbine3", ("3.367906049203724e-05", "[0.5626821330545027, 11.434005458108444]")), + ("kepler0", ("4.866123310307557e-05", "[20.802505677616598, 95.98217598131895]")), + ("kepler1", ("0.00021301499593723758", "[-248.33198742510353, -31.775362835450746]")), + ("kepler2", ("0.0012060524690583414", "[-645.8671575, 1288.34451875]")), + ("himmilbeau", ("0.0005369187050519044", "[-0.1416015625, 890.0]")) + ) + + val smtAA_Float64: Map[String, (String, String)] = Map( ("doppler", ("4.1911988101104756e-13", "[-137.64316836703838, -0.02944244059231351]")), - ("sine", ("7.410949021433787e-16", "[-1.0043283753002443, 1.0043283753002443]")), - ("sineOrder3", ("1.3400508287924788e-15", "[-1.0055351041384715, 1.0055351041384715]")), - ("sqroot", ("3.208544541166703e-13", "[-334.7934013614431, 1.6171710207127035]")), + ("sine", ("6.951079086011496e-16", "[-1.0043283753002443, 1.0043283753002443]")), + ("sineOrder3", ("1.2254703612493457e-15", "[-1.0055351041384715, 1.0055351041384715]")), + ("sqroot", ("3.088640454507186e-13", "[-334.7934013614431, 1.6171710207127035]")), ("bspline0", ("1.6190752442450204e-16", "[0.0, 0.16666666666666666]")), - ("bspline1", ("8.141635513917815e-16", "[0.16617838541666666, 0.6671549479166666]")), - ("bspline2", ("7.262708952756233e-16", "[0.16666666666666666, 0.6668650309244791]")), + ("bspline1", ("7.031412489292659e-16", "[0.16617838541666666, 0.6671549479166666]")), + ("bspline2", ("6.152485928131076e-16", "[0.16666666666666666, 0.6668650309244791]")), ("bspline3", ("1.0639637319324417e-16", "[-0.16666666666666666, 0.0]")), ("rigidBody1", ("3.2152058793144533e-13", "[-705.0, 705.0]")), ("rigidBody2", ("3.652100843964945e-11", "[-56049.1667175293, 58740.0]")), - ("verhulst", ("5.741159583498406e-16", "[0.3640144188500088, 0.9473239405662036]")), - ("predatorPrey", ("2.511128651395663e-16", "[0.03727705922396188, 0.33711264367110555]")), + ("verhulst", ("5.404475598160268e-16", "[0.3640144188500088, 0.9473239405662036]")), + ("predatorPrey", ("2.4756633989832674e-16", "[0.03727705922396188, 0.33711264367110555]")), ("carbonGas", ("3.7244002681234605e-08", "[4301713.35625, 1.6740286809375E7]")), - ("turbine1", ("9.140869202455534e-14", "[-18.536291471169722, -1.986569412987437]")), - ("turbine2", ("1.285668665546453e-13", "[-28.564652726384942, 3.8296552477598733]")), - ("turbine3", ("6.545577699383315e-14", "[0.5626821330545027, 11.434005458108444]")) + ("turbine1", ("8.868501650477874e-14", "[-18.536291471169722, -1.986569412987437]")), + ("turbine2", ("1.231195155150921e-13", "[-28.564652726384942, 3.8296552477598733]")), + ("turbine3", ("6.273210147405655e-14", "[0.5626821330545027, 11.434005458108444]")), + ("kepler0", ("9.063860773039779e-14", "[20.802505677616598, 95.98217598131895]")), + ("kepler1", ("4.251930363352586e-13", "[-278.40975267291554, -31.96861380296662]")), + ("kepler2", ("2.2464476501227186e-12", "[-645.8671575, 1288.34451875]")), + ("himmilbeau", ("1.0000889005823412e-12", "[-0.1416015625, 890.0]")) ) - test("interval-affine") { + + + + test("warm-up") { val ctx = context.copy(options=Seq(ChoiceOption("rangeMethod", "interval"))) val (_, program) = pipeline.run(ctx, inputPrg.deepCopy) + ctx.reporter.info("time taken (for all benchmarks): " + ctx.timers.analysis) + + } + + ignore("interval-interval") { + //TODO + } + + test("interval-affine-fixed32") { + val ctx = context.copy(options=Seq(ChoiceOption("rangeMethod", "interval"), + ChoiceOption("precision", "Fixed32"))) + val (resCtx, program) = pipeline.run(ctx, inputPrg.deepCopy) + ctx.reporter.info("time taken (for all benchmarks): " + ctx.timers.analysis) for (fnc <- program.defs) { - getLastExpression(fnc.body.get) match { - case x: NumAnnotation if x.hasError => - val absError = x.absError - val range = x.interval - assert(intAA(fnc.id.toString)._1 === absError.toString, fnc.id) - assert(intAA(fnc.id.toString)._2 === range.toString, fnc.id) - //println(s"${fnc.id} abs error: ${x.absError}, range: ${x.interval}") - - // if there is no error, then this is a bug - case x: NumAnnotation => - assert(false) - } + val absError = resCtx.resultAbsoluteErrors(fnc.id) + val range = resCtx.resultRealRanges(fnc.id) + + assert(intAA_Fixed32(fnc.id.toString)._1 === absError.toString, fnc.id) + assert(intAA_Fixed32(fnc.id.toString)._2 === range.toString, fnc.id) + } + } + + test("interval-affine-float32") { + val ctx = context.copy(options=Seq(ChoiceOption("rangeMethod", "interval"), + ChoiceOption("precision", "Float32"))) + val (resCtx, program) = pipeline.run(ctx, inputPrg.deepCopy) + ctx.reporter.info("time taken (for all benchmarks): " + ctx.timers.analysis) + + for (fnc <- program.defs) { + val absError = resCtx.resultAbsoluteErrors(fnc.id) + val range = resCtx.resultRealRanges(fnc.id) + assert(intAA_Float32(fnc.id.toString)._1 === absError.toString, fnc.id) + assert(intAA_Float32(fnc.id.toString)._2 === range.toString, fnc.id) + } + } + + test("interval-affine-float64") { + val ctx = context.copy(options=Seq(ChoiceOption("rangeMethod", "interval"), + ChoiceOption("precision", "Float64"))) + val (resCtx, program) = pipeline.run(ctx, inputPrg.deepCopy) + ctx.reporter.info("time taken (for all benchmarks): " + ctx.timers.analysis) + + for (fnc <- program.defs) { + val absError = resCtx.resultAbsoluteErrors(fnc.id) + val range = resCtx.resultRealRanges(fnc.id) + + assert(intAA_Float64(fnc.id.toString)._1 === absError.toString, fnc.id) + assert(intAA_Float64(fnc.id.toString)._2 === range.toString, fnc.id) } } - test("affine-affine") { - val ctx = context.copy(options=Seq(ChoiceOption("rangeMethod", "affine"))) - val (_, program) = pipeline.run(ctx, inputPrg.deepCopy) + test("affine-affine-fixed32") { + val ctx = context.copy(options=Seq(ChoiceOption("rangeMethod", "affine"), + ChoiceOption("precision", "Fixed32"))) + val (resCtx, program) = pipeline.run(ctx, inputPrg.deepCopy) + ctx.reporter.info("time taken (for all benchmarks): " + ctx.timers.analysis) + + for (fnc <- program.defs) { + val absError = resCtx.resultAbsoluteErrors(fnc.id) + val range = resCtx.resultRealRanges(fnc.id) + + assert(aaAA_Fixed32(fnc.id.toString)._1 === absError.toString, fnc.id) + assert(aaAA_Fixed32(fnc.id.toString)._2 === range.toString, fnc.id) + } + } + + test("affine-affine-float32") { + val ctx = context.copy(options=Seq(ChoiceOption("rangeMethod", "affine"), + ChoiceOption("precision", "Float32"))) + val (resCtx, program) = pipeline.run(ctx, inputPrg.deepCopy) + ctx.reporter.info("time taken (for all benchmarks): " + ctx.timers.analysis) for (fnc <- program.defs) { - getLastExpression(fnc.body.get) match { - case x: NumAnnotation if x.hasError => - val absError = x.absError - val range = x.interval - assert(aaAA(fnc.id.toString)._1 === absError.toString, fnc.id) - assert(aaAA(fnc.id.toString)._2 === range.toString, fnc.id) - //println(s"${fnc.id} abs error: ${x.absError}, range: ${x.interval}") - - // if there is no error, then this is a bug - case x: NumAnnotation => - assert(false) + val absError = resCtx.resultAbsoluteErrors(fnc.id) + val range = resCtx.resultRealRanges(fnc.id) + + assert(aaAA_Float32(fnc.id.toString)._1 === absError.toString, fnc.id) + assert(aaAA_Float32(fnc.id.toString)._2 === range.toString, fnc.id) + } + } + + test("affine-affine-float64") { + val ctx = context.copy(options=Seq(ChoiceOption("rangeMethod", "affine"), + ChoiceOption("precision", "Float64"))) + val (resCtx, program) = pipeline.run(ctx, inputPrg.deepCopy) + ctx.reporter.info("time taken (for all benchmarks): " + ctx.timers.analysis) + + for (fnc <- program.defs) { + val absError = resCtx.resultAbsoluteErrors(fnc.id) + val range = resCtx.resultRealRanges(fnc.id) + assert(aaAA_Float64(fnc.id.toString)._1 === absError.toString, fnc.id) + assert(aaAA_Float64(fnc.id.toString)._2 === range.toString, fnc.id) + } + } + + test("smt-affine-fixed32") { + val ctx = context.copy(options=Seq(ChoiceOption("rangeMethod", "smt"), + ChoiceOption("precision", "Fixed32"))) + val (resCtx, program) = pipeline.run(ctx, inputPrg.deepCopy) + ctx.reporter.info("time taken (for all benchmarks): " + ctx.timers.analysis) + + for (fnc <- program.defs) { + // kepler1 sometimes returns a different result, presumably because of a Z3 timeout + if (fnc.id.toString != "kepler1") { + val absError = resCtx.resultAbsoluteErrors(fnc.id) + val range = resCtx.resultRealRanges(fnc.id) + + assert(smtAA_Fixed32(fnc.id.toString)._1 === absError.toString, fnc.id) + assert(smtAA_Fixed32(fnc.id.toString)._2 === range.toString, fnc.id) } } } - // TODO: there seems to be an issue with the classpath, - // we get a NullpointerException at Solver.scala, line 37 - ignore("smt-affine") { - val ctx = context.copy(options=Seq(ChoiceOption("rangeMethod", "smt"))) - val (_, program) = pipeline.run(ctx, inputPrg.deepCopy) + test("smt-affine-float32") { + val ctx = context.copy(options=Seq(ChoiceOption("rangeMethod", "smt"), + ChoiceOption("precision", "Float32"))) + val (resCtx, program) = pipeline.run(ctx, inputPrg.deepCopy) + ctx.reporter.info("time taken (for all benchmarks): " + ctx.timers.analysis) for (fnc <- program.defs) { - getLastExpression(fnc.body.get) match { - case x: NumAnnotation if x.hasError => - val absError = x.absError - val range = x.interval - //assert(aaAA(fnc.id.toString)._1 === absError.toString) - //assert(aaAA(fnc.id.toString)._2 === range.toString) - //println(s"${fnc.id} abs error: ${x.absError}, range: ${x.interval}") - - // if there is no error, then this is a bug - case x: NumAnnotation => + // kepler1 sometimes returns a different result, presumably because of a Z3 timeout + if (fnc.id.toString != "kepler1") { + val absError = resCtx.resultAbsoluteErrors(fnc.id) + val range = resCtx.resultRealRanges(fnc.id) + assert(smtAA_Float32(fnc.id.toString)._1 === absError.toString, fnc.id) + assert(smtAA_Float32(fnc.id.toString)._2 === range.toString, fnc.id) + } + } + } + + + test("smt-affine-float64") { + val ctx = context.copy(options=Seq(ChoiceOption("rangeMethod", "smt"), + ChoiceOption("precision", "Float64"))) + val (resCtx, program) = pipeline.run(ctx, inputPrg.deepCopy) + ctx.reporter.info("time taken (for all benchmarks): " + ctx.timers.analysis) + + for (fnc <- program.defs) { + // kepler1 sometimes returns a different result, presumably because of a Z3 timeout + if (fnc.id.toString != "kepler1") { + val absError = resCtx.resultAbsoluteErrors(fnc.id) + val range = resCtx.resultRealRanges(fnc.id) + assert(smtAA_Float64(fnc.id.toString)._1 === absError.toString, fnc.id) + assert(smtAA_Float64(fnc.id.toString)._2 === range.toString, fnc.id) } } } @@ -159,47 +422,31 @@ class AbsErrorAnalysisRegressionTest extends FunSuite { { val ctx = jetContext.copy(options=Seq(ChoiceOption("rangeMethod", "affine"))) - val (_, prg) = pipeline.run(ctx, inputPrgJet.deepCopy) + val (resCtx, prg) = pipeline.run(ctx, inputPrgJet.deepCopy) for (fnc <- prg.defs) { - getLastExpression(fnc.body.get) match { - case x: NumAnnotation if x.hasError => - val absError = x.absError - val range = x.interval - assert("4.088949761228701e-08" === absError.toString) - assert("[-2236030.8006656803, 2243586.0576923075]" === range.toString) - //println(s"${fnc.id} abs error: ${x.absError}, range: ${x.interval}") - - // if there is no error, then this is a bug - case x: NumAnnotation => - assert(false) - } + val absError = resCtx.resultAbsoluteErrors(fnc.id) + val range = resCtx.resultRealRanges(fnc.id) + assert("4.088949761228701e-08" === absError.toString) + assert("[-2236030.8006656803, 2243586.0576923075]" === range.toString) } } { val ctx = jetContext.copy(options=Seq(ChoiceOption("rangeMethod", "smt"))) - val (_, prg) = pipeline.run(ctx, inputPrgJet.deepCopy) + val (resCtx, prg) = pipeline.run(ctx, inputPrgJet.deepCopy) for (fnc <- prg.defs) { - getLastExpression(fnc.body.get) match { - case x: NumAnnotation if x.hasError => - val absError = x.absError - val range = x.interval - assert("1.1589127503449215e-08" === absError.toString) - assert("[-1781.3056274865728, 4820.245415085443]" === range.toString) - //println(s"${fnc.id} abs error: ${x.absError}, range: ${x.interval}") - - // if there is no error, then this is a bug - case x: NumAnnotation => - assert(false) - } + val absError = resCtx.resultAbsoluteErrors(fnc.id) + val range = resCtx.resultRealRanges(fnc.id) + assert("1.1589127503449215e-08" === absError.toString) + assert("[-1781.3056274865728, 4820.245415085443]" === range.toString) } } } - private def getLastExpression(e: Expr): Expr = e match { - case Let(_, _, body) => getLastExpression(body) - case _ => e - } + // private def getLastExpression(e: Expr): Expr = e match { + // case Let(_, _, body) => getLastExpression(body) + // case _ => e + // } } diff --git a/src/test/scala/regression/FixedpointCodegenRegressionTest.scala b/src/test/scala/regression/FixedpointCodegenRegressionTest.scala new file mode 100644 index 0000000000000000000000000000000000000000..1eb55a18dbd4514b7e86ba7afd8f27297b487232 --- /dev/null +++ b/src/test/scala/regression/FixedpointCodegenRegressionTest.scala @@ -0,0 +1,225 @@ + +package daisy +package test.regression + +import scala.collection.immutable.Seq +import org.scalatest.FunSuite + +import utils.FinitePrecision._ + +/** + Regression test for the fixed-point arithmetic code generation. +*/ +class FixedpointCodegenRegressionTest extends FunSuite { + + val inputFile: String = "src/test/resources/AbsErrorRegressionFunctions.scala" + + val context = Context( + reporter = new DefaultReporter(Set()), + options = List(), + files = List(inputFile) + ) + + val inputPrg = frontend.ExtractionPhase(context) + val pipeline = analysis.SpecsProcessingPhase andThen + transform.SSATransformerPhase andThen + analysis.RangeErrorPhase + + test("codegen-doppler-fixed32") { + val ctx = context.copy(options=Seq( + ChoiceOption("rangeMethod", "interval"), + ChoiceOption("precision", "Fixed32"), + ListOption("functions", List("doppler")))) + val (_, program) = pipeline.run(ctx, inputPrg.deepCopy) + + val fnc = program.defs.find(_.id.toString == "doppler").get + //println(fnc.body.get) + val newBody = backend.CodeGenerationPhase.toFixedPointCode(fnc.body.get, Fixed(32)) + //println(newBody) + + assert(newBody.toString === + "(let (_tmp := ((2576980378 * T) >> 31)) in\n" + + " (let (t1 := (((2779984691 << 4) + _tmp) >> 4)) in\n" + + " (let (_tmp1 := -(t1)) in\n" + + " (let (_tmp4 := ((_tmp1 * v) >> 31)) in\n" + + " (let (_tmp2 := (((t1 << 2) + u) >> 2)) in\n" + + " (let (_tmp3 := (((t1 << 2) + u) >> 2)) in\n" + + " (let (_tmp5 := ((_tmp2 * _tmp3) >> 32)) in\n" + + " ((_tmp4 << 29) / _tmp5))))))))" + ) + } + + test("codegen-sine-fixed32") { + val ctx = context.copy(options=Seq( + ChoiceOption("rangeMethod", "interval"), + ChoiceOption("precision", "Fixed32"), + ListOption("functions", List("sine")))) + val (_, program) = pipeline.run(ctx, inputPrg.deepCopy) + + val fnc = program.defs.find(_.id.toString == "sine").get + //println(fnc.body.get) + + val newBody = backend.CodeGenerationPhase.toFixedPointCode(fnc.body.get, Fixed(32)) + + assert(newBody.toString === + "(let (_tmp229 := ((x * x) >> 32)) in\n" + + " (let (_tmp230 := ((_tmp229 * x) >> 31)) in\n" + + " (let (_tmp231 := ((_tmp230 << 31) / 3221225472)) in\n" + + " (let (_tmp236 := (((x << 1) - _tmp231) >> 2)) in\n" + + " (let (_tmp232 := ((x * x) >> 32)) in\n" + + " (let (_tmp233 := ((_tmp232 * x) >> 31)) in\n" + + " (let (_tmp234 := ((_tmp233 * x) >> 32)) in\n" + + " (let (_tmp235 := ((_tmp234 * x) >> 32)) in\n" + + " (let (_tmp237 := ((_tmp235 << 29) / 4026531840)) in\n" + + " (let (_tmp244 := (((_tmp236 << 2) + _tmp237) >> 2)) in\n" + + " (let (_tmp238 := ((x * x) >> 32)) in\n" + + " (let (_tmp239 := ((_tmp238 * x) >> 31)) in\n" + + " (let (_tmp240 := ((_tmp239 * x) >> 32)) in\n" + + " (let (_tmp241 := ((_tmp240 * x) >> 32)) in\n" + + " (let (_tmp242 := ((_tmp241 * x) >> 31)) in\n" + + " (let (_tmp243 := ((_tmp242 * x) >> 32)) in\n" + + " (let (_tmp245 := ((_tmp243 << 24) / 2642411520)) in\n" + + " (((_tmp244 << 2) - _tmp245) >> 2))))))))))))))))))" + ) + } + + test("codegen-rigidBody1-fixed32") { + val ctx = context.copy(options=Seq( + ChoiceOption("rangeMethod", "interval"), + ChoiceOption("precision", "Fixed32"), + ListOption("functions", List("rigidBody1")))) + val (_, program) = pipeline.run(ctx, inputPrg.deepCopy) + + val fnc = program.defs.find(_.id.toString == "rigidBody1").get + + val newBody = backend.CodeGenerationPhase.toFixedPointCode(fnc.body.get, Fixed(32)) + //println(newBody) + + assert(newBody.toString === + "(let (_tmp510 := -(x1)) in\n" + + " (let (_tmp512 := ((_tmp510 * x2) >> 32)) in\n" + + " (let (_tmp511 := ((2147483648 * x2) >> 31)) in\n" + + " (let (_tmp513 := ((_tmp511 * x3) >> 32)) in\n" + + " (let (_tmp514 := ((_tmp512 - (_tmp513 << 1)) >> 2)) in\n" + + " (let (_tmp515 := (((_tmp514 << 6) - x1) >> 6)) in\n" + + " (((_tmp515 << 6) - x3) >> 6)))))))") + } + + test("codegen-turbine1-fixed32") { + val ctx = context.copy(options=Seq( + ChoiceOption("rangeMethod", "interval"), + ChoiceOption("precision", "Fixed32"), + ListOption("functions", List("turbine1")))) + val (_, program) = pipeline.run(ctx, inputPrg.deepCopy) + + val fnc = program.defs.find(_.id.toString == "turbine1").get + + val newBody = backend.CodeGenerationPhase.toFixedPointCode(fnc.body.get, Fixed(32)) + //println(newBody) + + assert(newBody.toString === + "(let (_tmp771 := ((r * r) >> 32)) in\n" + + " (let (_tmp772 := ((2147483648 << 28) / _tmp771)) in\n" + + " (let (_tmp781 := (((3221225472 << 2) + _tmp772) >> 2)) in\n" + + " (let (_tmp773 := ((2147483648 * v) >> 31)) in\n" + + " (let (_tmp774 := ((3221225472 - (_tmp773 << 2)) >> 2)) in\n" + + " (let (_tmp777 := ((536870912 * _tmp774) >> 29)) in\n" + + " (let (_tmp775 := ((w * w) >> 32)) in\n" + + " (let (_tmp776 := ((_tmp775 * r) >> 32)) in\n" + + " (let (_tmp778 := ((_tmp776 * r) >> 32)) in\n" + + " (let (_tmp779 := ((_tmp777 * _tmp778) >> 32)) in\n" + + " (let (_tmp780 := ((2147483648 - (v << 2)) >> 2)) in\n" + + " (let (_tmp782 := ((_tmp779 << 30) / _tmp780)) in\n" + + " (let (_tmp783 := ((_tmp781 - (_tmp782 << 4)) >> 4)) in\n" + + " (((_tmp783 << 3) - 2415919104) >> 3))))))))))))))") + } + + test("codegen-turbine2-fixed32") { + val ctx = context.copy(options=Seq( + ChoiceOption("rangeMethod", "interval"), + ChoiceOption("precision", "Fixed32"), + ListOption("functions", List("turbine2")))) + val (_, program) = pipeline.run(ctx, inputPrg.deepCopy) + + val fnc = program.defs.find(_.id.toString == "turbine2").get + + val newBody = backend.CodeGenerationPhase.toFixedPointCode(fnc.body.get, Fixed(32)) + //println(newBody) + + assert(newBody.toString === + "(let (_tmp1013 := ((3221225472 * v) >> 31)) in\n" + + " (let (_tmp1009 := ((2147483648 * v) >> 31)) in\n" + + " (let (_tmp1007 := ((w * w) >> 32)) in\n" + + " (let (_tmp1008 := ((_tmp1007 * r) >> 32)) in\n" + + " (let (_tmp1010 := ((_tmp1008 * r) >> 32)) in\n" + + " (let (_tmp1011 := ((_tmp1009 * _tmp1010) >> 31)) in\n" + + " (let (_tmp1012 := ((2147483648 - (v << 2)) >> 2)) in\n" + + " (let (_tmp1014 := ((_tmp1011 << 29) / _tmp1012)) in\n" + + " (let (_tmp1015 := ((_tmp1013 - (_tmp1014 << 2)) >> 2)) in\n" + + " (((_tmp1015 << 5) - 2684354560) >> 5))))))))))") + + } + + + test("codegen-kepler0-fixed32") { + val ctx = context.copy(options=Seq( + ChoiceOption("rangeMethod", "interval"), + ChoiceOption("precision", "Fixed32"), + ListOption("functions", List("kepler0")))) + val (_, program) = pipeline.run(ctx, inputPrg.deepCopy) + + val fnc = program.defs.find(_.id.toString == "kepler0").get + + backend.CodeGenerationPhase.reporter = ctx.reporter + val newBody = backend.CodeGenerationPhase.toFixedPointCode(fnc.body.get, Fixed(32)) + //println(newBody) + + assert(newBody.toString === + "(let (_tmp1252 := ((x2 * x5) >> 32)) in\n" + + " (let (_tmp1253 := ((x3 * x6) >> 32)) in\n" + + " (let (_tmp1254 := ((_tmp1252 + _tmp1253) >> 1)) in\n" + + " (let (_tmp1255 := ((x2 * x3) >> 32)) in\n" + + " (let (_tmp1256 := (((_tmp1254 << 1) - _tmp1255) >> 1)) in\n" + + " (let (_tmp1257 := ((x5 * x6) >> 32)) in\n" + + " (let (_tmp1264 := ((_tmp1256 << 1) - _tmp1257)) in\n" + + " (let (_tmp1258 := -(x1)) in\n" + + " (let (_tmp1259 := ((_tmp1258 + x2) << 1)) in\n" + + " (let (_tmp1260 := ((_tmp1259 + (x3 << 1)) >> 2)) in\n" + + " (let (_tmp1261 := ((_tmp1260 << 1) - x4)) in\n" + + " (let (_tmp1262 := ((_tmp1261 + x5) >> 1)) in\n" + + " (let (_tmp1263 := (((_tmp1262 << 1) + x6) >> 2)) in\n" + + " (let (_tmp1265 := ((x1 * _tmp1263) >> 31)) in\n" + + " ((_tmp1264 + (_tmp1265 << 1)) >> 2)))))))))))))))") + } + + test("codegen-himmilbeau-fixed32") { + val ctx = context.copy(options=Seq( + ChoiceOption("rangeMethod", "interval"), + ChoiceOption("precision", "Fixed32"), + ListOption("functions", List("himmilbeau")))) + val (_, program) = pipeline.run(ctx, inputPrg.deepCopy) + + val fnc = program.defs.find(_.id.toString == "himmilbeau").get + + val newBody = backend.CodeGenerationPhase.toFixedPointCode(fnc.body.get, Fixed(32)) + //println(newBody) + + assert(newBody.toString === + "(let (_tmp1547 := ((x1 * x1) >> 31)) in\n" + + " (let (_tmp1548 := (((_tmp1547 << 2) + x2) >> 2)) in\n" + + " (let (_tmp1551 := (((_tmp1548 << 1) - 2952790016) >> 2)) in\n" + + " (let (_tmp1549 := ((x1 * x1) >> 31)) in\n" + + " (let (_tmp1550 := (((_tmp1549 << 2) + x2) >> 2)) in\n" + + " (let (_tmp1552 := (((_tmp1550 << 1) - 2952790016) >> 2)) in\n" + + " (let (_tmp1559 := ((_tmp1551 * _tmp1552) >> 31)) in\n" + + " (let (_tmp1553 := ((x2 * x2) >> 31)) in\n" + + " (let (_tmp1554 := ((x1 + (_tmp1553 << 2)) >> 2)) in\n" + + " (let (_tmp1557 := (((_tmp1554 << 2) - 3758096384) >> 3)) in\n" + + " (let (_tmp1555 := ((x2 * x2) >> 31)) in\n" + + " (let (_tmp1556 := ((x1 + (_tmp1555 << 2)) >> 2)) in\n" + + " (let (_tmp1558 := (((_tmp1556 << 2) - 3758096384) >> 3)) in\n" + + " (let (_tmp1560 := ((_tmp1557 * _tmp1558) >> 31)) in\n" + + " ((_tmp1559 + _tmp1560) >> 1)))))))))))))))") + } + +}