Commit c8a22163 authored by Heiko Becker's avatar Heiko Becker
Browse files

Fix error bounds for simpleDoppler example, correct parameter bound...

Fix error bounds for simpleDoppler example, correct parameter bound computations and implemented addition bounds
parent 159cf56b
Require Import Coq.QArith.QArith Coq.QArith.Qminmax Coq.QArith.Qabs Coq.ZArith.ZArith Coq.Reals.Reals Coq.QArith.Qreals Coq.QArith.Qcanon.
Require Import Coq.QArith.QArith Coq.QArith.Qminmax Coq.QArith.Qabs Coq.ZArith.ZArith Coq.Reals.Reals Coq.QArith.Qreals.
(*Coq.QArith.Qcanon.*)
Require Import Coq.micromega.Psatz.
Require Import Daisy.Expressions Daisy.Infra.RationalConstruction Daisy.Infra.RealSimps.
Require Import Daisy.Expressions Daisy.Infra.RationalConstruction Daisy.Infra.RealSimps Daisy.IntervalArith.
Section ComputableErrors.
......@@ -8,6 +9,7 @@ Section ComputableErrors.
Definition error :Type := Q.
Definition analysisResult :Type := exp Q -> interval * error.
(*
Definition Qc2Q (q:Qc) :Q := match q with
Qcmake a P => a end.
......@@ -16,7 +18,7 @@ Section ComputableErrors.
Proof.
intros q. unfold Q2Qc.
unfold Qc2Q; auto.
Qed.
Qed. *)
Definition Qleb :=Qle_bool.
......@@ -27,6 +29,9 @@ Section ComputableErrors.
Definition maxAbs (iv:interval) :=
Qmax (Qabs (fst iv)) (Qabs (snd iv)).
Definition RmaxAbsFun (iv:interval) :=
Rmax (Rabs (Q2R (fst iv))) (Rabs (Q2R (snd iv))).
Lemma maxAbs_pointInterval a:
maxAbs (a,a) == Qabs a.
Proof.
......@@ -65,21 +70,29 @@ Section ComputableErrors.
Definition machineEpsilon := (1#(2^53)).
Fixpoint validErrorbound (e:exp Q) (env:analysisResult) :=
match e with
|Var v => false
|Param v => false
(* A constant will be a point interval. Take maxAbs never the less *)
|Const n => let (intv,err) := env (Const n) in
Qleb (maxAbs intv * machineEpsilon) err
|Binop b e1 e2 =>
let rec := andb (validErrorbound e1 env) (validErrorbound e2 env) in
match b with
Plus => false
|Sub => false
|Mult => false
|Div => false
end
end.
let (intv, err) := (env e) in
match e with
|Var v => false
|Param v => Qleb (maxAbs intv * machineEpsilon) err
(* A constant will be a point interval. Take maxAbs never the less *)
|Const n => Qleb (maxAbs intv * machineEpsilon) err
|Binop b e1 e2 =>
let (ive1, err1) := env e1 in
let (ive2, err2) := env e2 in
let rec := andb (validErrorbound e1 env) (validErrorbound e2 env) in
let upperBoundE1 := maxAbs ive1 in
let upperBoundE2 := maxAbs ive2 in
let e1F := upperBoundE1 + (upperBoundE1 * machineEpsilon) in
let e2F := upperBoundE2 + (upperBoundE2 * machineEpsilon) in
let theVal :=
match b with
| Plus => Qleb (err1 + err2 + (Qabs (e1F + e2F) * machineEpsilon)) err
| Sub => false
| Mult => false
| Div => false
end
in andb rec theVal
end.
Definition u:nat := 1.
(** 1655/5 = 331; 0,4 = 2/5 **)
......@@ -91,11 +104,11 @@ Section ComputableErrors.
Definition valCstAddVarU:exp Q := Binop Plus valCst varU.
(** Error values **)
Definition errCst1 :Q := rationalFromNum (7358558207215537#1) 15 14.
Definition errVaru := rationalFromNum (2220446049250313#1) 15 14.
Definition lowerBoundAddUCst:R := 1157 / 5.
Definition upperBoundAddUCst:R := 2157 / 5.
Definition errAddUCst := rationalFromNum (19158008512931703#1) 16 13.
Definition errCst1 :Q := rationalFromNum (36792791036077685#1) 16 14.
Definition errVaru := rationalFromNum (11102230246251565#1) 16 14.
Definition lowerBoundAddUCst:Q := 1157 # 5.
Definition upperBoundAddUCst:Q := 2157 # 5.
Definition errAddUCst := rationalFromNum (9579004256465851#1) 15 14.
(** The added assertion becomes the precondition for us **)
Definition precondition := fun env:nat -> R => (- 100 <= env u)%R /\ (env u <= 100)%R.
......@@ -106,7 +119,25 @@ Section ComputableErrors.
|theVar: absEnv varU (mkInterval (- 100) (100)) errVaru
|theAddition: absEnv valCstAddVarU (mkInterval lowerBoundAddUCst upperBoundAddUCst) errAddUCst. **)
Definition validConstant := Eval compute in validErrorbound (valCst) (fun x => ((cst1,cst1),errCst1)).
Definition validConstant := Eval compute in validErrorbound (valCst) (fun x => ((cst1,cst1),errCst1)).
Definition validParam := Eval compute in validErrorbound (varU) (fun x => (((-100) # 1,100 # 1),errVaru)).
Definition validAddition := Eval compute in validErrorbound (valCstAddVarU)
(fun x => match x with
|Param _ => ((cst1,cst1),errCst1)
|Const _ => (((-100) # 1,100 # 1),errVaru)
|_ => ((lowerBoundAddUCst,upperBoundAddUCst),errAddUCst) end).
Definition envApproximatesPrecond (P:nat -> interval) (absenv:analysisResult) :=
forall u:nat,
let (ivlo,ivhi) := P u in
let (absIv,err) := absenv (Param Q u) in
let (abslo,abshi) := absIv in
(abslo <= ivlo /\ ivhi <= abshi)%Q.
Definition precondValidForExec (P:nat->interval) (cenv:nat->R) :=
forall v:nat,
let (ivlo,ivhi) := P v in
(Q2R ivlo <= cenv v <= Q2R ivhi)%R.
Lemma Q2R0_is_0:
Q2R 0 = 0%R.
......@@ -115,7 +146,7 @@ Section ComputableErrors.
Qed.
Lemma Rabs_eq_Qabs:
forall x, Rabs (Q2R x) = Q2R (Qabs x).
forall x, Rabs (Q2R x) = Q2R (Qabs x).
Proof.
intro.
apply Qabs_case; unfold Rabs; destruct Rcase_abs; intros; try auto.
......@@ -144,6 +175,41 @@ Section ComputableErrors.
rewrite Ropp_0; auto.
Qed.
Lemma maxAbs_impl_RmaxAbs (iv:interval):
RmaxAbsFun iv = Q2R (maxAbs iv).
Proof.
unfold RmaxAbsFun, maxAbs.
repeat rewrite Rabs_eq_Qabs.
apply Q.max_case_strong.
- intros x y x_eq_y Rmax_x.
rewrite Rmax_x.
unfold Q2R. rewrite <- RMicromega.Rinv_elim.
setoid_rewrite Rmult_comm at 1.
+ rewrite <- Rmult_assoc.
rewrite <- RMicromega.Rinv_elim.
rewrite <- mult_IZR.
rewrite x_eq_y.
rewrite mult_IZR.
rewrite Rmult_comm; auto.
+ simpl.
hnf; intros.
pose proof (pos_INR_nat_of_P (Qden y)).
rewrite H in H0.
lra.
+ simpl; hnf; intros.
pose proof (pos_INR_nat_of_P (Qden x)).
rewrite H in H0; lra.
- intros snd_le_fst.
apply Qle_Rle in snd_le_fst.
apply Rmax_left in snd_le_fst.
subst; auto.
- intros fst_le_snd.
apply Qle_Rle in fst_le_snd.
apply Rmax_right in fst_le_snd.
subst; auto.
Qed.
Lemma validErrorboundCorrectConstant:
forall cenv absenv (n:Q) nR nF e,
eval_exp 0%R cenv (Const (Q2R n)) nR ->
......@@ -175,6 +241,55 @@ Section ComputableErrors.
apply error_valid.
Qed.
Lemma validErrorboundCorrectParam:
forall cenv absenv (v:nat) nR nF e P plo phi,
envApproximatesPrecond P absenv ->
precondValidForExec P cenv ->
eval_exp 0%R cenv (Param R v) nR ->
eval_exp (Q2R machineEpsilon) cenv (Param R v) nF ->
validErrorbound (Param Q v) absenv = true ->
absenv (Param Q v) = ((plo,phi),e) ->
(Rabs (nR - nF) <= (Q2R e))%R.
Proof.
intros cenv absenv v nR nF e P plo phi absenv_approx_p cenv_approx_p eval_real eval_float error_valid absenv_param.
unfold validErrorbound in error_valid.
rewrite absenv_param in error_valid.
inversion eval_real; subst.
rewrite perturb_0_val in eval_real; auto.
rewrite perturb_0_val; auto.
inversion eval_float; subst.
unfold envApproximatesPrecond in absenv_approx_p.
unfold precondValidForExec in cenv_approx_p.
specialize (absenv_approx_p v).
specialize (cenv_approx_p v).
assert (exists ivlo ivhi, (ivlo,ivhi) = P v) by (destruct (P v) as [ivlo ivhi]; repeat eexists).
destruct H as [ivlo [ivhi P_eq]].
rewrite <- P_eq in absenv_approx_p, cenv_approx_p.
rewrite absenv_param in absenv_approx_p.
destruct absenv_approx_p as [plo_le_ivlo ivhi_le_phi].
destruct cenv_approx_p as [ivlo_le_cenv cenv_le_ivhi].
apply Qle_Rle in plo_le_ivlo; apply Qle_Rle in ivhi_le_phi.
pose proof (Rle_trans (Q2R plo) (Q2R ivlo) (cenv v) plo_le_ivlo ivlo_le_cenv).
pose proof (Rle_trans (cenv v) (Q2R ivhi) (Q2R phi) cenv_le_ivhi ivhi_le_phi).
pose proof (RmaxAbs (Q2R plo) (cenv v) (Q2R phi) H H2).
unfold perturb.
rewrite Rabs_err_simpl.
rewrite Rabs_mult.
eapply Rle_trans.
eapply Rmult_le_compat; [ apply Rabs_pos | apply Rabs_pos | | ].
(* TODO: Does that work out? *)
apply H3.
apply H1.
pose proof (maxAbs_impl_RmaxAbs (plo,phi)).
unfold RmaxAbsFun in H4.
simpl in H4.
rewrite H4.
apply Qle_bool_iff in error_valid.
apply Qle_Rle in error_valid.
rewrite Q2R_mult in error_valid.
apply error_valid.
Qed.
Lemma validErrorboundCorrect:
forall cenv (n:Q),
eval_exp 0%R cenv (Const n%R) nR ->
......
......@@ -4,6 +4,29 @@ Require Import Daisy.AbsoluteError Daisy.Commands Daisy.IntervalArith
Daisy.Expressions Daisy.ErrorBounds
Daisy.Infra.RealConstruction Daisy.Infra.Abbrevs Daisy.Infra.RealSimps.
(*
TODO: update according to:
[ Info ]
[ Info ] Starting specs preprocessing phase
[ Info ] Finished specs preprocessing phase
[ Info ]
[ Info ]
[ Info ] Starting range-error phase
[ Info ] Machine epsilon 1.1102230246251565E-16
[ Info ] 331.4: [331.4, 331.4],3.6792791036077685e-14
[ Info ] u: [-100.0, 100.0],1.1102230246251565e-14
[ Info ] (331.4 + u): [231.4, 431.4],9.579004256465851e-14
[ Info ] Finished range-error phase
[ Info ]
[ Info ] Starting info phase
[ Info ] doppler
[ Info ] abs-error: 9.579004256465851e-14, range: [231.4, 431.4],
[ Info ] rel-error: 4.139586973407887E-16
[ Info ] Finished info phase
[ Info ] time:
[ Info ] info: 7 ms, rangeError: 67 ms, analysis: 11 ms, frontend: 2171 ms,
[ Info ]
[ Info ] Starting specs preprocessing phase
[ Info ] Finished specs preprocessing phase
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment