Commit 176780f1 authored by Heiko Becker's avatar Heiko Becker

Try to prove copnstant case, fix machine epsilon. Still not working...

parent bf0162f3
......@@ -32,7 +32,7 @@ Inductive exp (V:Type): Type :=
Define the machine epsilon for floating point operations.
FIXME: Currently set to 1.0 instead of the concrete value!
**)
Definition machineEpsilon:R := realFromNum 59604645 7 8.
Definition machineEpsilon:R := realFromNum 2220446049250313 15 16.
(**
Define a perturbation function to ease writing of basic definitions
**)
......
Require Import Coq.Reals.Reals Interval.Interval_tactic.
Require Import Coq.Reals.Reals Interval.Interval_tactic Coq.micromega.Psatz.
Require Import Daisy.abs_err Daisy.daisy_lang Daisy.newIntervalArith Daisy.exps Daisy.realConversion.
(*
Daisy output:
[ Info ]
[ Info ] Starting specs preprocessing phase
[ Info ] Finished specs preprocessing phase
[ Info ]
[ Info ]
[ Info ] Starting range-error phase
[ Info ] Machine epsilon 5.9604645E-8
[ Info ] 331.4: [331.4, 331.4], 2.8421709430404007e-14
[ Info ] u: [-100.0, 100.0], 7.105427357601002e-15
[ Info ] (331.4 + u): [231.4, 431.4],6.394884621840902e-14
[ Info ] Finished range-error phase
[ Info ]
[ Info ] Starting info phase
[ Info ] doppler
[ Info ] abs-error: 6.394884621840902e-14, range: [231.4, 431.4],
[ Info ] rel-error: 2.7635629307869066E-16
[ Info ] Finished info phase
[ Info ] time:
[ Info ] info: 8 ms, rangeError: 81 ms, analysis: 14 ms, frontend: 2547 ms,
def doppler(u: Real): Real = {
require(-100.0 <= u && u <= 100)
331.4 + u
}
[ Info ]
[ Info ] Starting specs preprocessing phase
[ Info ] Finished specs preprocessing phase
[ Info ]
[ Info ]
[ Info ] Starting range-error phase
[ Info ] Machine epsilon 2.220446049250313E-16
[ Info ] 331.4: [331.4, 331.4],7.358558207215537e-14
[ Info ] u: [-100.0, 100.0],2.220446049250313e-14
[ Info ] (331.4 + u): [231.4, 431.4],1.9158008512931703e-13
[ Info ] Finished range-error phase
[ Info ]
[ Info ] Starting info phase
[ Info ] doppler
[ Info ] abs-error: 1.9158008512931703e-13, range: [231.4, 431.4],
[ Info ] rel-error: 8.279173946815775E-16
[ Info ] Finished info phase
[ Info ] time:
[ Info ] info: 7 ms, rangeError: 58 ms, analysis: 8 ms, frontend: 2167 ms,
*)
Definition u:nat := 1.
(** 1655/5 = 331; 0,4 = 2/5 **)
......@@ -41,11 +32,11 @@ Definition valCst:exp R := Const cst1.
Definition valCstAddVarU:exp R := Binop Plus valCst varU.
(** Error values **)
Definition errCst1 := realFromNum 28421709430404007 16 14.
Definition errVaru := realFromNum 7105427357601002 15 15.
Definition errCst1 := realFromNum 7358558207215537 15 14.
Definition errVaru := realFromNum 2220446049250313 15 14.
Definition lowerBoundAddUCst:R := - (1157) / 5.
Definition upperBoundAddUCst:R := 2157 / 5.
Definition errAddUCst := realFromNum 6394884621840902 15 14.
Definition errAddUCst := realFromNum 19158008512931703 16 13.
(** The added assertion becomes the precondition for us **)
Definition precondition := fun env:nat -> R => (- 100 <= env u)%R /\ (env u <= 100)%R.
......@@ -62,7 +53,7 @@ Inductive absEnv : abs_env :=
Goal forall (cenv:nat->R) (vR:R) (vF:R),
precondition cenv ->
eval_exp (0)%R cenv valCstAddVarU vR ->
eval_exp m_eps cenv valCstAddVarU vF ->
eval_exp machineEpsilon cenv valCstAddVarU vF ->
AbsErrExp valCstAddVarU absEnv errAddUCst /\ (Rabs (vR - vF) <= errAddUCst)%R.
Proof.
intros cenv vR vF precond eval_real eval_float.
......@@ -76,10 +67,35 @@ Proof.
+ unfold valCst.
apply (AbsErrConst cst1 (mkInterval cst1 cst1) errCst1); [constructor | ].
unfold isSoundErr; simpl.
admit. (** Should be easy once we have the real epsilon **)
unfold errCst1, cst1, machineEpsilon.
unfold realFromNum.
unfold negPow.
unfold Rabs.
assert (1657 / 5 - 2220446049250313 * (1 / 10 ^ 15) * (1 / 2 ^ 16) >= 0)%R by interval.
assert (1657 / 5 + 2220446049250313 * (1 / 10 ^ 15) * (1 / 2 ^ 16) >= 0)%R by interval.
destruct Rcase_abs.
* exfalso.
apply Rlt_not_le in r.
apply r.
apply Rge_le in H; auto.
* destruct Rcase_abs.
{ exfalso.
apply Rlt_not_le in r0.
apply r0.
apply Rge_le in H0; auto.
}
{
assert ((1657 / 5 - 2220446049250313 * (1 / 10 ^ 15) * (1 / 2 ^ 16)) <= 1657 / 5 + 2220446049250313 * (1 / 10 ^ 15) * (1 / 2 ^ 16))%R
by interval.
assert (Rmax (1657 / 5 - 2220446049250313 * (1 / 10 ^ 15) * (1 / 2 ^ 16))
(1657 / 5 + 2220446049250313 * (1 / 10 ^ 15) * (1 / 2 ^ 16)) = (1657 / 5 + 2220446049250313 * (1 / 10 ^ 15) * (1 / 2 ^ 16)))%R by (apply Rmax_right; auto).
rewrite H2.
field_simplify.
(** TODO: Flaw in computation. Wolfram Alpha says that this inequation is simply false **)
admit.
+ apply (AbsErrVar u (mkInterval (- 100) (100)) errVaru); [constructor | ].
unfold isSoundErr; simpl.
unfold m_eps, errVaru.
unfold machineEpsilon, errVaru.
unfold realFromNum.
unfold negPow.
simpl.
......@@ -131,4 +147,4 @@ Proof.
as floatAddErrContained
by (apply additionIsValid; auto).
admit.
Admitted.
\ No newline at end of file
Admitted.
......@@ -40,7 +40,9 @@ object RangeErrorPhase extends DaisyPhase with ErrorFunctions {
override def run(ctx: Context, prg: Program): (Context, Program) = {
reporter = ctx.reporter
reporter.info(s"\nStarting $name")
reporter.info(s"Machine epsilon ${math.ulp(1.0f)/2}")
val eps = java.lang.Math.ulp(1.0)
reporter.info(s"Machine epsilon $eps")
val timer = ctx.timers.rangeError.start
......
......@@ -42,7 +42,8 @@ object FinitePrecision {
}
def absRoundoff(r: Rational): Rational = {
// PERFORMANCE: this may not be the fastest way
Rational.fromDouble(math.ulp(Rational.abs(r).floatValue)/2)
Rational.fromDouble(math.ulp(1.0))*Rational.abs(r)
//Rational.fromDouble(math.ulp(Rational.abs(r).floatValue)/2)
}
}
......@@ -58,7 +59,8 @@ object FinitePrecision {
}
def absRoundoff(r: Rational): Rational = {
Rational.fromDouble(math.ulp(Rational.abs(r).doubleValue)/2)
Rational.fromDouble(math.ulp(1.0))*Rational.abs(r)
//Rational.fromDouble(math.ulp(Rational.abs(r).doubleValue)/2)
}
}
......@@ -103,4 +105,4 @@ object FinitePrecision {
val minNormal: Rational = zero // dummy
}*/
}
\ No newline at end of file
}
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