mirror of
https://github.com/AdaCore/why3.git
synced 2026-02-12 12:34:55 -08:00
1087 lines
38 KiB
Plaintext
1087 lines
38 KiB
Plaintext
(** {1 Unbounded floating-point numbers}
|
|
|
|
|
|
See also {h <a href="https://inria.hal.science/hal-04343157">this report</a>}
|
|
|
|
*)
|
|
|
|
module UFloat
|
|
use real.RealInfix
|
|
use real.FromInt
|
|
use real.Abs
|
|
use ieee_float.RoundingMode
|
|
|
|
type t
|
|
|
|
val function uround mode real : t
|
|
val function to_real t : real
|
|
val function of_int int : t
|
|
axiom to_real_of_int : forall x [of_int x]. to_real (of_int x) = from_int x
|
|
|
|
val function uzero : t
|
|
axiom uzero_spec : to_real uzero = 0.0
|
|
|
|
val function uone : t
|
|
axiom uone_spec : to_real uone = 1.0
|
|
|
|
val function utwo : t
|
|
axiom utwo_spec : to_real utwo = 2.0
|
|
|
|
constant eps:real
|
|
constant eta:real
|
|
axiom eps_bounds : 0. <. eps <. 1.
|
|
axiom eta_bounds : 0. <. eta <. 1.
|
|
|
|
(* To avoid "inline_trivial" to break the forward_propagation strategy *)
|
|
meta "inline:no" function eps
|
|
meta "inline:no" function eta
|
|
|
|
let function uadd (x y:t) : t
|
|
(* TODO: Do we want the two first assertions in our context ?
|
|
We only use them to prove the addition lemma *)
|
|
ensures { abs (to_real result -. (to_real x +. to_real y)) <=. abs (to_real x) }
|
|
ensures { abs (to_real result -. (to_real x +. to_real y)) <=. abs (to_real y) }
|
|
ensures {
|
|
abs (to_real result -. (to_real x +. to_real y))
|
|
<=. abs (to_real x +. to_real y) *. eps
|
|
}
|
|
= uround RNE (to_real x +. to_real y)
|
|
|
|
let function usub (x y:t) : t
|
|
(* TODO: Do we want the two first assertions in our context ?
|
|
We only use them to prove the addition lemma *)
|
|
ensures { abs (to_real result -. (to_real x -. to_real y)) <=. abs (to_real x) }
|
|
ensures { abs (to_real result -. (to_real x -. to_real y)) <=. abs (to_real y) }
|
|
ensures {
|
|
abs (to_real result -. (to_real x -. to_real y))
|
|
<=. abs (to_real x -. to_real y) *. eps
|
|
}
|
|
= uround RNE (to_real x -. to_real y)
|
|
|
|
let function umul (x y:t) : t
|
|
ensures {
|
|
abs (to_real result -. (to_real x *. to_real y))
|
|
<=. abs (to_real x *. to_real y) *. eps +. eta
|
|
}
|
|
= uround RNE (to_real x *. to_real y)
|
|
|
|
let function udiv (x y:t) : t
|
|
requires { to_real y <> 0. }
|
|
ensures {
|
|
abs (to_real result -. (to_real x /. to_real y))
|
|
<=. abs (to_real x /. to_real y) *. eps +. eta
|
|
}
|
|
= uround RNE (to_real x /. to_real y)
|
|
|
|
let function uminus (x:t) : t
|
|
ensures { to_real result = -. (to_real x) }
|
|
= uround RNE (-. (to_real x))
|
|
|
|
predicate is_exact (uop : t -> t -> t) (x y :t)
|
|
|
|
(* Exact division but can still underflow, giving eta as error *)
|
|
let function udiv_exact (x y:t) : t
|
|
requires { to_real y <> 0. }
|
|
requires { is_exact udiv x y }
|
|
ensures { abs (to_real result -. (to_real x /. to_real y)) <=. eta }
|
|
= uround RNE (to_real x /. to_real y)
|
|
|
|
(** Infix operators *)
|
|
let function ( ++. ) (x:t) (y:t) : t = uadd x y
|
|
let function ( --. ) (x:t) (y:t) : t = usub x y
|
|
let function ( **. ) (x:t) (y:t) : t = umul x y
|
|
(* Why3 doesn't support abbreviations so we need to add the requires *)
|
|
let function ( //. ) (x:t) (y:t) : t
|
|
requires { to_real y <> 0. }
|
|
= udiv x y
|
|
let function ( --._ ) (x:t) : t = uminus x
|
|
let function ( ///. ) (x:t) (y:t) : t
|
|
requires { to_real y <> 0. }
|
|
requires { is_exact udiv x y }
|
|
= udiv_exact x y
|
|
|
|
(* Some constants *)
|
|
constant u0:t
|
|
axiom to_real_u0 : to_real u0 = 0.0
|
|
constant u1:t
|
|
axiom to_real_u1 : to_real u1 = 1.0
|
|
constant u2:t
|
|
axiom to_real_u2 : to_real u2 = 2.0
|
|
constant u4:t
|
|
axiom to_real_u4 : to_real u4 = 4.0
|
|
constant u8:t
|
|
axiom to_real_u8 : to_real u8 = 8.0
|
|
constant u16:t
|
|
axiom to_real_u16 : to_real u16 = 16.0
|
|
constant u32:t
|
|
axiom to_real_u32 : to_real u32 = 32.0
|
|
constant u64:t
|
|
axiom to_real_u64 : to_real u64 = 64.0
|
|
constant u128:t
|
|
axiom to_real_u128 : to_real u128 = 128.0
|
|
constant u256:t
|
|
axiom to_real_u256 : to_real u256 = 256.0
|
|
constant u512:t
|
|
axiom to_real_u512 : to_real u512 = 512.0
|
|
constant u1024:t
|
|
axiom to_real_u1024 : to_real u1024 = 1024.0
|
|
constant u2048:t
|
|
axiom to_real_u2048 : to_real u2048 = 2048.0
|
|
constant u4096:t
|
|
axiom to_real_u4096 : to_real u4096 = 4096.0
|
|
constant u8192:t
|
|
axiom to_real_u8192 : to_real u8192 = 8192.0
|
|
constant u16384:t
|
|
axiom to_real_u16384 : to_real u16384 = 16384.0
|
|
constant u32768:t
|
|
axiom to_real_u32768 : to_real u32768 = 32768.0
|
|
constant u65536:t
|
|
axiom to_real_u65536 : to_real u65536 = 65536.0
|
|
|
|
predicate is_positive_power_of_2 (x:t) =
|
|
x = u1 \/ x = u2 || x = u4 || x = u8 || x = u16 || x = u32 || x = u64
|
|
|| x = u128 \/ x = u256 || x = u4096 || x = u8192 || x = u16384 || x = u32768
|
|
|| x = u65536
|
|
|
|
axiom div_by_positive_power_of_2_is_exact :
|
|
forall x y. is_positive_power_of_2 y -> is_exact udiv x y
|
|
end
|
|
|
|
(** {2 Single-precision unbounded floats} *)
|
|
module USingle
|
|
use real.RealInfix
|
|
|
|
type usingle
|
|
|
|
constant eps:real = 0x1p-24 /. (1. +. 0x1p-24)
|
|
constant eta:real = 0x1p-150
|
|
|
|
clone export UFloat with
|
|
type t = usingle,
|
|
constant eps = eps,
|
|
constant eta = eta,
|
|
axiom.
|
|
end
|
|
|
|
|
|
(** {3 Double-precision unbounded floats} *)
|
|
module UDouble
|
|
use real.RealInfix
|
|
type udouble
|
|
|
|
constant eps:real = 0x1p-53 /. (1. +. 0x1p-53)
|
|
constant eta:real = 0x1p-1075
|
|
|
|
clone export UFloat with
|
|
type t = udouble,
|
|
constant eps = eps,
|
|
constant eta = eta,
|
|
axiom.
|
|
end
|
|
|
|
(* Helper lemmas to help the proof of propagation lemmas *)
|
|
module HelperLemmas
|
|
use real.RealInfix
|
|
use real.Abs
|
|
|
|
let ghost div_order_compat (x y z:real)
|
|
requires { x <=. y }
|
|
requires { 0. <. z }
|
|
ensures { x /. z <=. y /. z }
|
|
= ()
|
|
|
|
let ghost div_order_compat2 (x y z:real)
|
|
requires { x <=. y }
|
|
requires { 0. >. z }
|
|
ensures { y /. z <=. x /. z }
|
|
= ()
|
|
|
|
let ghost mult_err (x x_exact x_factor x_rel x_abs y:real)
|
|
requires { 0. <=. x_rel }
|
|
requires { 0. <=. x_abs }
|
|
requires { abs x_exact <=. x_factor }
|
|
requires { abs (x -. x_exact) <=. x_rel *. x_factor +. x_abs }
|
|
ensures { abs (x *. y -. x_exact *. y) <=. x_rel *. abs (x_factor *. y) +. x_abs *. abs y }
|
|
=
|
|
assert {
|
|
y >=. 0. ->
|
|
abs (x *. y -. x_exact *. y) <=. abs (x_rel *. x_factor *. y) +. x_abs *. abs y
|
|
by
|
|
(x_exact -. x_rel *. x_factor -. x_abs) *. y <=. x *. y <=. (x_exact +. x_rel *. x_factor +. x_abs) *. y
|
|
};
|
|
assert {
|
|
y <. 0. ->
|
|
abs (x *. y -. x_exact *. y) <=. abs (x_rel *. x_factor *. y) +. x_abs *. abs y
|
|
by
|
|
(x_exact +. x_rel *. x_factor +. x_abs) *. y <=. x *. y <=. (x_exact -. x_rel *. x_factor -. x_abs) *. y
|
|
}
|
|
|
|
let ghost mult_err_combine (x x_exact x_factor x_rel x_abs y exact_y y_factor y_rel y_abs:real)
|
|
requires { 0. <=. x_rel }
|
|
requires { 0. <=. y_rel }
|
|
requires { 0. <=. x_abs }
|
|
requires { 0. <=. y_abs }
|
|
requires { abs x_exact <=. x_factor }
|
|
requires { abs exact_y <=. y_factor }
|
|
requires { abs (x -. x_exact) <=. x_rel *. x_factor +. x_abs }
|
|
requires { abs (y -. exact_y) <=. y_rel *. y_factor +. y_abs }
|
|
ensures {
|
|
abs (x *. y -. x_exact *. exact_y)
|
|
<=. (x_rel +. y_rel +. x_rel *. y_rel) *. (x_factor *. y_factor)
|
|
+. (y_abs +. y_abs *. x_rel) *. x_factor
|
|
+. (x_abs +. x_abs *. y_rel) *. y_factor
|
|
+. x_abs *. y_abs
|
|
}
|
|
=
|
|
mult_err x x_exact x_factor x_rel x_abs y;
|
|
mult_err y exact_y y_factor y_rel y_abs x_exact;
|
|
mult_err y exact_y y_factor y_rel y_abs x_factor;
|
|
assert {
|
|
abs (x *. y -. x_exact *. exact_y) <=. (y_rel *. x_factor *. y_factor) +. (y_abs *. x_factor) +. (x_rel *. abs (x_factor *. y)) +. x_abs *. abs y
|
|
};
|
|
assert {
|
|
abs (x *. y -. x_exact *. exact_y) <=. (y_rel *. x_factor *. y_factor) +. (x_rel *. (x_factor *. y_factor *. (1. +. y_rel) +. x_factor *. y_abs)) +. y_abs *. x_factor +. x_abs *. abs y
|
|
by
|
|
abs (x_factor *. y) <=. x_factor *. y_factor *. (1. +. y_rel) +. x_factor *. y_abs
|
|
};
|
|
assert {
|
|
x_abs *. abs y <=. x_abs *. (y_factor *. (1. +. y_rel) +. y_abs)
|
|
}
|
|
|
|
use real.ExpLog
|
|
|
|
let ghost exp_approx_err (x x_approx x_factor a b :real)
|
|
requires { abs (x_approx -. x) <=. x_factor *. a +. b }
|
|
requires { x <=. x_factor }
|
|
ensures {
|
|
abs (exp(x_approx) -. exp(x)) <=. exp(x) *. (exp(a *. x_factor +. b) -. 1.)
|
|
}
|
|
=
|
|
assert {
|
|
exp(x_approx) <=. exp(x) +. exp(x) *. (exp(a *. x_factor +. b) -. 1.)
|
|
by
|
|
exp (x_approx) <=. exp(x) *. exp (a *. x_factor +. b)
|
|
};
|
|
assert {
|
|
exp(x_approx) >=. exp(x) -. exp(x) *. (exp(a *. x_factor +. b) -. 1.)
|
|
by
|
|
exp (x_approx) >=. exp(x) *. exp (-. a *. x_factor -. b)
|
|
so
|
|
exp(x_approx) -. exp(x) >=. exp(x) *. (exp (-. a *. x_factor -. b) -. 1.)
|
|
so
|
|
exp(a *. x_factor +. b) +. exp(-.a *. x_factor -. b) >=. 2.
|
|
so
|
|
-. exp(a *. x_factor +. b) +. 1. <=. exp(-.a *. x_factor -. b) -. 1.
|
|
so
|
|
exp(x) *. ((-. exp(a *. x_factor +. b)) +. 1.) <=. exp(x) *. (exp(-. a *. x_factor -. b) -. 1.)
|
|
so
|
|
-. exp(x) *. (exp(a *. x_factor +. b) -. 1.) <=. exp(x) *. (exp(-. a *. x_factor -. b) -. 1.)
|
|
}
|
|
|
|
let lemma log_1_minus_x (x:real)
|
|
requires { 0. <=. abs x <. 1. }
|
|
ensures { log (1. +. x) <=. -. log (1. -. x) }
|
|
=
|
|
assert { 1. +. x <=. 1. /. (1. -. x) };
|
|
assert { forall x y z. 0. <=. x -> 0. <. y -> 0. <=. z -> x *. y <=. z -> x <=. z /. y };
|
|
assert { exp (-.log (1. -. x)) = 1. /. (1. -. x) }
|
|
|
|
let lemma log2_1_minus_x (x:real)
|
|
requires { 0. <=. abs x <. 1. }
|
|
ensures { log2 (1. +. x) <=. -. log2 (1. -. x) }
|
|
=
|
|
div_order_compat (log (1. +. x)) (-. log (1. -. x)) (log 2.);
|
|
log_1_minus_x x
|
|
|
|
let lemma log10_1_minus_x (x:real)
|
|
requires { 0. <=. abs x <. 1. }
|
|
ensures { log10 (1. +. x) <=. -. log10 (1. -. x) }
|
|
=
|
|
div_order_compat (log (1. +. x)) (-. log (1. -. x)) (log 10.);
|
|
log_1_minus_x x
|
|
|
|
let ghost log_approx_err (x x_approx x_factor a b :real)
|
|
requires { abs (x_approx -. x) <=. x_factor *. a +. b }
|
|
requires { 0. <. (x -. a *. x_factor -. b) }
|
|
requires { 0. <. x <=. x_factor }
|
|
ensures {
|
|
abs (log x_approx -. log x) <=. -. log(1. -. ((a *. x_factor +. b) /. x))
|
|
}
|
|
=
|
|
assert { a *. x_factor +. b = x *. ((a *. x_factor +. b) /. x) };
|
|
assert {
|
|
log (x *. (1. -. (a *. x_factor +. b) /. x))
|
|
<=. log x_approx
|
|
<=. log (x *. (1. +. (a *. x_factor +. b) /.x))
|
|
by
|
|
0. <. (x -. (a *. x_factor +. b)) <=. x_approx
|
|
};
|
|
log_1_minus_x ((a *. x_factor +. b) /. x)
|
|
|
|
let ghost log2_approx_err (x x_approx x_factor a b :real)
|
|
requires { abs (x_approx -. x) <=. x_factor *. a +. b }
|
|
requires { 0. <. (x -. a *. x_factor -. b) }
|
|
requires { 0. <. x <=. x_factor }
|
|
ensures {
|
|
abs (log2 x_approx -. log2 x) <=. -. log2(1. -. ((a *. x_factor +. b) /. x))
|
|
}
|
|
=
|
|
assert { a *. x_factor +. b = x *. ((a *. x_factor +. b) /. x) };
|
|
assert {
|
|
log2 (x *. (1. -. (a *. x_factor +. b) /. x))
|
|
<=. log2 x_approx
|
|
<=. log2 (x *. (1. +. (a *. x_factor +. b) /.x))
|
|
by
|
|
0. <. (x -. (a *. x_factor +. b)) <=. x_approx
|
|
};
|
|
log2_1_minus_x ((a *. x_factor +. b) /. x)
|
|
|
|
let ghost log10_approx_err (x x_approx x_factor a b :real)
|
|
requires { abs (x_approx -. x) <=. x_factor *. a +. b }
|
|
requires { 0. <. (x -. a *. x_factor -. b) }
|
|
requires { 0. <. x <=. x_factor }
|
|
ensures {
|
|
abs (log10 x_approx -. log10 x) <=. -. log10(1. -. ((a *. x_factor +. b) /. x))
|
|
}
|
|
=
|
|
assert { a *. x_factor +. b = x *. ((a *. x_factor +. b) /. x) };
|
|
assert {
|
|
log10 (x *. (1. -. (a *. x_factor +. b) /. x))
|
|
<=. log10 x_approx
|
|
<=. log10 (x *. (1. +. (a *. x_factor +. b) /.x))
|
|
by
|
|
0. <. (x -. (a *. x_factor +. b)) <=. x_approx
|
|
};
|
|
log10_1_minus_x ((a *. x_factor +. b) /. x)
|
|
|
|
use real.Trigonometry
|
|
|
|
lemma sin_of_approx : forall x y. abs (sin x -. sin y) <=. abs (x -. y)
|
|
lemma cos_of_approx : forall x y. abs (cos x -. cos y) <=. abs (x -. y)
|
|
|
|
use real.Sum
|
|
use int.Int
|
|
use real.FromInt
|
|
|
|
let rec ghost sum_approx_err (fi_rel fi_abs:real) (f f_exact f_factor : int -> real) (a b:int)
|
|
requires { a <= b }
|
|
requires { forall i. a <= i < b -> abs (f i -. f_exact i) <=. f_factor i *. fi_rel +. fi_abs }
|
|
variant { b - a }
|
|
ensures { abs (sum f a b -. sum f_exact a b) <=. fi_rel *. sum f_factor a b +. fi_abs *. from_int (b-a) }
|
|
=
|
|
if (a < b) then
|
|
begin
|
|
sum_approx_err fi_rel fi_abs f f_exact f_factor a (b - 1)
|
|
end
|
|
|
|
end
|
|
|
|
(** {4 Single propagation lemmas} *)
|
|
module USingleLemmas
|
|
use real.RealInfix
|
|
use real.FromInt
|
|
use real.Abs
|
|
use USingle
|
|
|
|
let lemma uadd_single_error_propagation (x_f y_f r: usingle) (x x_factor x_rel x_abs y y_factor y_rel y_abs : real)
|
|
requires {
|
|
abs (to_real x_f -. x) <=. x_rel *. x_factor +. x_abs
|
|
}
|
|
requires {
|
|
abs (to_real y_f -. y) <=. y_rel *. y_factor +. y_abs
|
|
}
|
|
requires { abs x <=. x_factor }
|
|
requires { abs y <=. y_factor }
|
|
(* TODO: Use (0 <=. x_rel \/ (x_factor = 0 /\ x_abs = 0)), same for y. *)
|
|
requires { 0. <=. x_rel }
|
|
requires { 0. <=. y_rel }
|
|
requires { 0. <=. x_abs }
|
|
requires { 0. <=. y_abs }
|
|
requires { r = (x_f ++. y_f) }
|
|
ensures {
|
|
abs (to_real r -. (x +. y)) <=.
|
|
(x_rel +. y_rel +. eps) *. (x_factor +. y_factor)
|
|
+. ((1. +. eps +. y_rel) *. x_abs +. (1. +. eps +. x_rel) *. y_abs)
|
|
}
|
|
=
|
|
let ghost delta = abs (to_real (x_f ++. y_f) -. (to_real x_f +. to_real y_f)) in
|
|
assert {
|
|
0. <=. x_rel /\ 0. <=. y_rel ->
|
|
delta <=.
|
|
(eps +. y_rel) *. x_factor +. (eps +. x_rel) *. y_factor
|
|
+. (x_rel +. eps) *. y_abs +. (y_rel +. eps) *. x_abs
|
|
by
|
|
(delta <=. x_factor *. x_rel +. x_abs +. x_factor
|
|
so
|
|
x_factor +. x_abs <=. eps *. (y_factor +. y_abs) ->
|
|
delta <=. (eps +. x_rel) *. y_factor
|
|
+. (eps +. y_rel) *. x_factor
|
|
+. (y_rel +. eps) *. x_abs +. (x_rel +. eps) *. y_abs
|
|
by
|
|
delta <=. eps *. (y_factor +. y_abs) *. x_rel
|
|
+. (eps *. (y_factor +. y_abs)))
|
|
/\
|
|
(delta <=. y_factor *. y_rel +. y_abs +. y_factor
|
|
so
|
|
abs y_factor +. y_abs <=. eps *. (x_factor +. x_abs) ->
|
|
delta <=. (eps +. y_rel) *. x_factor
|
|
+. (eps +. x_rel) *. y_factor
|
|
+. (x_rel +. eps) *. y_abs +. (y_rel +. eps) *. x_abs
|
|
by
|
|
delta <=. eps *. (x_factor +. x_abs) *. y_rel
|
|
+. (eps *. (x_factor +. x_abs)))
|
|
/\
|
|
(
|
|
(eps *. (x_factor +. x_abs) <. abs y_factor +. y_abs /\
|
|
eps *. (y_factor +. y_abs) <. abs x_factor +. x_abs) ->
|
|
(delta <=.
|
|
(eps +. y_rel) *. x_factor +. (eps +. x_rel) *. y_factor
|
|
+. (x_rel +. eps) *. y_abs +. (y_rel +. eps) *. x_abs
|
|
by
|
|
abs (to_real x_f +. to_real y_f) <=.
|
|
abs (to_real x_f -. x) +. x_factor +. (abs (to_real y_f -. y) +. y_factor)
|
|
so
|
|
x_factor *. x_rel <=. (y_factor +. y_abs) /. eps *. x_rel /\
|
|
y_factor *. y_rel <=. (x_factor +. x_abs) /. eps *. y_rel))
|
|
}
|
|
|
|
let lemma usub_single_error_propagation (x_f y_f r : usingle) (x x_factor x_rel x_abs y y_factor y_rel y_abs : real)
|
|
requires {
|
|
abs (to_real x_f -. x) <=. x_rel *. x_factor +. x_abs
|
|
}
|
|
requires {
|
|
abs (to_real y_f -. y) <=. y_rel *. y_factor +. y_abs
|
|
}
|
|
requires { abs x <=. x_factor }
|
|
requires { abs y <=. y_factor }
|
|
requires { 0. <=. x_abs }
|
|
requires { 0. <=. y_abs }
|
|
requires { 0. <=. x_rel }
|
|
requires { 0. <=. y_rel }
|
|
requires { r = x_f --. y_f }
|
|
ensures {
|
|
abs (to_real r -. (x -. y))
|
|
<=. (x_rel +. y_rel +. eps) *. (x_factor +. y_factor)
|
|
+. ((1. +. eps +. y_rel) *. x_abs +. (1. +. eps +. x_rel) *. y_abs)
|
|
}
|
|
= uadd_single_error_propagation x_f (--. y_f) r x x_factor x_rel x_abs (-. y) y_factor y_rel y_abs
|
|
|
|
use HelperLemmas
|
|
|
|
let lemma umul_single_error_propagation (x_f y_f r : usingle) (x x_factor x_rel x_abs y y_factor y_rel y_abs : real)
|
|
requires {
|
|
abs (to_real x_f -. x) <=. x_rel *. x_factor +. x_abs
|
|
}
|
|
requires {
|
|
abs (to_real y_f -. y) <=. y_rel *. y_factor +. y_abs
|
|
}
|
|
requires { abs x <=. x_factor }
|
|
requires { abs y <=. y_factor }
|
|
requires { 0. <=. x_rel }
|
|
requires { 0. <=. y_rel }
|
|
requires { 0. <=. x_abs }
|
|
requires { 0. <=. y_abs }
|
|
requires { r = x_f **. y_f }
|
|
ensures {
|
|
abs (to_real r -. (x *. y)) <=.
|
|
(eps +. (x_rel +. y_rel +. x_rel *. y_rel) *. (1. +. eps)) *. (x_factor *. y_factor)
|
|
+. (((y_abs +. y_abs *. x_rel) *. x_factor
|
|
+. (x_abs +. x_abs *. y_rel) *. y_factor
|
|
+. x_abs *. y_abs) *. (1. +. eps) +. eta)
|
|
}
|
|
=
|
|
assert {
|
|
to_real x_f *. to_real y_f -. abs (to_real x_f *. to_real y_f) *. eps -. eta
|
|
<=. to_real (x_f **. y_f)
|
|
<=. to_real x_f *. to_real y_f +. abs (to_real x_f *. to_real y_f) *. eps +. eta
|
|
};
|
|
assert { abs (x *. y) <=. x_factor *. y_factor by
|
|
abs x *. abs y <=. x_factor *. abs y = abs y *. x_factor <=. y_factor *. x_factor };
|
|
mult_err_combine (to_real x_f) x x_factor x_rel x_abs (to_real y_f) y y_factor y_rel y_abs
|
|
|
|
use real.ExpLog
|
|
|
|
let lemma log_single_error_propagation (logx_f x_f : usingle)
|
|
(x_exact x_factor log_rel log_abs x_rel x_abs : real)
|
|
requires { abs (to_real x_f -. x_exact) <=. x_rel *. x_factor +. x_abs }
|
|
requires {
|
|
abs (to_real logx_f -. log(to_real x_f))
|
|
<=. log_rel *. abs (log (to_real x_f)) +. log_abs
|
|
}
|
|
requires { 0. <. x_exact <=. x_factor }
|
|
requires { 0. <. (x_exact -. x_rel *. x_factor -. x_abs) }
|
|
requires { 0. <=. log_rel }
|
|
ensures {
|
|
abs (to_real logx_f -. log (x_exact))
|
|
<=. log_rel *. abs (log (x_exact)) +.
|
|
(-. log (1. -. ((x_rel *. x_factor +. x_abs) /. x_exact)) *. (1. +. log_rel)
|
|
+. log_abs)
|
|
}
|
|
=
|
|
log_approx_err x_exact (to_real x_f) x_factor x_rel x_abs;
|
|
assert {
|
|
abs (log (to_real x_f)) *. log_rel
|
|
<=. (abs (log (x_exact)) -. log (1.0 -. (((x_rel *. x_factor) +. x_abs) /. x_exact))) *. log_rel
|
|
}
|
|
|
|
let lemma log2_single_error_propagation (log2x_f x_f : usingle)
|
|
(x_exact x_factor log_rel log_abs x_rel x_abs : real)
|
|
requires { abs (to_real x_f -. x_exact) <=. x_rel *. x_factor +. x_abs }
|
|
requires {
|
|
abs (to_real log2x_f -. log2(to_real x_f))
|
|
<=. log_rel *. abs (log2 (to_real x_f)) +. log_abs
|
|
}
|
|
requires { 0. <. x_exact <=. x_factor }
|
|
requires { 0. <. (x_exact -. x_rel *. x_factor -. x_abs) }
|
|
requires { 0. <=. log_rel }
|
|
ensures {
|
|
abs (to_real log2x_f -. log2 (x_exact))
|
|
<=. log_rel *. abs (log2 (x_exact)) +.
|
|
(-. log2 (1. -. ((x_rel *. x_factor +. x_abs) /. x_exact)) *. (1. +. log_rel)
|
|
+. log_abs)
|
|
}
|
|
=
|
|
log2_approx_err x_exact (to_real x_f) x_factor x_rel x_abs;
|
|
assert {
|
|
abs (log2 (to_real x_f)) *. log_rel
|
|
<=. (abs (log2 (x_exact)) -. log2 (1.0 -. (((x_rel *. x_factor) +. x_abs) /. x_exact))) *. log_rel
|
|
}
|
|
|
|
let lemma log10_single_error_propagation (log10x_f x_f : usingle)
|
|
(x_exact x_factor log_rel log_abs x_rel x_abs : real)
|
|
requires { abs (to_real x_f -. x_exact) <=. x_rel *. x_factor +. x_abs }
|
|
requires {
|
|
abs (to_real log10x_f -. log10(to_real x_f))
|
|
<=. log_rel *. abs (log10 (to_real x_f)) +. log_abs
|
|
}
|
|
requires { 0. <. x_exact <=. x_factor }
|
|
requires { 0. <. (x_exact -. x_rel *. x_factor -. x_abs) }
|
|
requires { 0. <=. log_rel }
|
|
ensures {
|
|
abs (to_real log10x_f -. log10 (x_exact))
|
|
<=. log_rel *. abs (log10 (x_exact)) +.
|
|
(-. log10 (1. -. ((x_rel *. x_factor +. x_abs) /. x_exact)) *. (1. +. log_rel)
|
|
+. log_abs)
|
|
}
|
|
=
|
|
log10_approx_err x_exact (to_real x_f) x_factor x_rel x_abs;
|
|
assert {
|
|
abs (log10 (to_real x_f)) *. log_rel
|
|
<=. (abs (log10 (x_exact)) -. log10 (1.0 -. (((x_rel *. x_factor) +. x_abs) /. x_exact))) *. log_rel
|
|
}
|
|
|
|
let lemma exp_single_error_propagation (expx_f x_f : usingle)
|
|
(x_exact x_factor exp_rel exp_abs x_rel x_abs : real)
|
|
requires { abs (to_real x_f -. x_exact) <=. x_rel *. x_factor +. x_abs }
|
|
requires {
|
|
abs (to_real expx_f -. exp(to_real x_f))
|
|
<=. exp_rel *. exp (to_real x_f) +. exp_abs
|
|
}
|
|
requires { x_exact <=. x_factor }
|
|
requires { 0. <=. exp_rel <=. 1. }
|
|
ensures {
|
|
abs (to_real expx_f -. exp (x_exact))
|
|
<=. (exp_rel +. (exp(x_rel *. x_factor +. x_abs) -. 1.) *. (1. +. exp_rel)) *. exp(x_exact)
|
|
+. exp_abs
|
|
}
|
|
=
|
|
exp_approx_err x_exact (to_real x_f) x_factor x_rel x_abs;
|
|
assert {
|
|
exp x_exact *. (1. -. exp_rel) -.
|
|
exp x_exact *. (exp (x_rel *. x_factor +. x_abs) -. 1.) *. (1. -. exp_rel)
|
|
-. exp_abs
|
|
<=. to_real expx_f
|
|
by
|
|
(exp x_exact -. exp x_exact *. (exp (x_rel *. x_factor +. x_abs) -. 1.))
|
|
*. (1. -. exp_rel) -. exp_abs
|
|
<=. exp (to_real x_f) *. (1. -. exp_rel) -. exp_abs
|
|
<=. to_real expx_f
|
|
};
|
|
assert {
|
|
to_real expx_f <=. (exp(x_exact) +. exp(x_exact)*.(exp(x_rel *. x_factor +. x_abs) -. 1.))*. (1. +. exp_rel) +. exp_abs
|
|
by
|
|
to_real expx_f <=. exp(to_real x_f) *. (1. +. exp_rel) +. exp_abs
|
|
};
|
|
|
|
|
|
use real.Trigonometry
|
|
|
|
let lemma sin_single_error_propagation (sinx_f x_f : usingle)
|
|
(x_exact x_factor sin_rel sin_abs x_rel x_abs : real)
|
|
requires { abs (to_real x_f -. x_exact) <=. x_rel *. x_factor +. x_abs }
|
|
requires {
|
|
abs (to_real sinx_f -. sin(to_real x_f))
|
|
<=. sin_rel *. abs (sin (to_real x_f)) +. sin_abs
|
|
}
|
|
requires { x_exact <=. x_factor }
|
|
requires { 0. <=. sin_rel }
|
|
ensures {
|
|
abs (to_real sinx_f -. sin (x_exact))
|
|
<=. sin_rel *. abs(sin(x_exact))
|
|
+. (((x_rel *. x_factor +. x_abs) *. (1. +. sin_rel)) +. sin_abs)
|
|
}
|
|
=
|
|
assert {
|
|
abs (sin (to_real x_f)) *. sin_rel
|
|
<=. (abs (sin x_exact) +. (x_rel *. x_factor +. x_abs)) *. sin_rel
|
|
}
|
|
|
|
let lemma cos_single_error_propagation (cosx_f x_f : usingle)
|
|
(x_exact x_factor cos_rel cos_abs x_rel x_abs : real)
|
|
requires { abs (to_real x_f -. x_exact) <=. x_rel *. x_factor +. x_abs }
|
|
requires {
|
|
abs (to_real cosx_f -. cos(to_real x_f))
|
|
<=. cos_rel *. abs (cos (to_real x_f)) +. cos_abs
|
|
}
|
|
requires { x_exact <=. x_factor }
|
|
requires { 0. <=. cos_rel }
|
|
ensures {
|
|
abs (to_real cosx_f -. cos (x_exact))
|
|
<=. cos_rel *. abs(cos(x_exact))
|
|
+. (((x_rel *. x_factor +. x_abs) *. (1. +. cos_rel)) +. cos_abs)
|
|
}
|
|
=
|
|
assert {
|
|
abs (cos (to_real x_f)) *. cos_rel
|
|
<=. (abs (cos x_exact) +. (x_rel *. x_factor +. x_abs)) *. cos_rel
|
|
}
|
|
|
|
use real.Sum
|
|
use int.Int
|
|
use real.FromInt
|
|
|
|
function real_fun (f:int -> usingle) : int -> real = fun i -> to_real (f i)
|
|
|
|
let lemma sum_single_error_propagation (x : usingle)
|
|
(f : int -> usingle) (f_exact f_factor f_factor' : int -> real) (n:int)
|
|
(sum_rel sum_abs f_rel f_abs : real)
|
|
requires {
|
|
forall i. 0 <= i < n ->
|
|
abs ((real_fun f) i -. f_exact i) <=. f_rel *. f_factor i +. f_abs
|
|
}
|
|
requires {
|
|
forall i. 0 <= i < n ->
|
|
f_factor i -. f_rel *. f_factor i -. f_abs <=. f_factor' i <=. f_factor i +. f_rel *. f_factor i +. f_abs
|
|
}
|
|
requires {
|
|
abs (to_real x -. (sum (real_fun f) 0 n))
|
|
<=. sum_rel *. (sum f_factor' 0 n) +. sum_abs
|
|
}
|
|
requires { 0. <=. sum_rel }
|
|
requires { 0 <= n }
|
|
ensures {
|
|
abs (to_real x -. sum f_exact 0 n)
|
|
<=. (f_rel +. (sum_rel *. (1. +. f_rel))) *. sum f_factor 0 n +.
|
|
((f_abs *. from_int n *.(1. +. sum_rel)) +. sum_abs)
|
|
}
|
|
=
|
|
sum_approx_err f_rel f_abs (real_fun f) f_exact f_factor 0 n;
|
|
sum_approx_err f_rel f_abs f_factor' f_factor f_factor 0 n;
|
|
assert {
|
|
sum_rel *. sum f_factor' 0 n <=.
|
|
sum_rel *. (sum f_factor 0 n +. ((f_rel *. sum f_factor 0 n) +. (f_abs *. from_int n)))
|
|
}
|
|
|
|
(* We don't have an error on y_f because in practice we won't have an exact division with an approximate divisor *)
|
|
let lemma udiv_exact_single_error_propagation (x_f y_f r: usingle) (x x_factor x_rel x_abs : real)
|
|
requires {
|
|
abs (to_real x_f -. x) <=. x_rel *. x_factor +. x_abs
|
|
}
|
|
requires { abs x <=. x_factor }
|
|
requires { 0. <=. x_rel }
|
|
requires { 0. <=. x_abs }
|
|
requires { 0. <> to_real y_f }
|
|
requires { is_exact udiv x_f y_f }
|
|
requires { r = x_f ///. y_f }
|
|
ensures {
|
|
abs (to_real r -. (x /. (to_real y_f))) <=.
|
|
x_rel *. (x_factor /. abs (to_real y_f)) +. ((x_abs /. abs (to_real y_f)) +. eta)
|
|
}
|
|
=
|
|
let lemma y_f_pos ()
|
|
requires { 0. <. to_real y_f }
|
|
ensures {
|
|
abs (to_real r -. (x /. (to_real y_f))) <=.
|
|
x_rel *. (x_factor /. to_real y_f) +. ((x_abs /. to_real y_f) +. eta)
|
|
}
|
|
=
|
|
div_order_compat (to_real x_f) (x +. x_rel *. x_factor +. x_abs) (to_real y_f);
|
|
div_order_compat (x -. x_rel *. x_factor -. x_abs) (to_real x_f) (to_real y_f) in
|
|
let lemma y_f_neg ()
|
|
requires { to_real y_f <. 0. }
|
|
ensures {
|
|
abs (to_real r -. (x /. (to_real y_f))) <=.
|
|
x_rel *. (x_factor /. abs (to_real y_f)) +. ((x_abs /. abs (to_real y_f)) +. eta)
|
|
}
|
|
=
|
|
div_order_compat2 (to_real x_f) (x +. x_rel *. x_factor +. x_abs) (to_real y_f);
|
|
(* TODO: Prove this somehow *)
|
|
assert {
|
|
forall x y. y <> 0.0 -> x /. y <=. abs x /. abs y
|
|
by abs (x /. y) = abs (x *. inv y) = abs x *. abs (inv y) = abs x *. inv (abs y) = abs x /. abs y
|
|
};
|
|
assert {
|
|
(x -. x_rel *. x_factor -. x_abs) /. to_real y_f
|
|
<=. x /. (to_real y_f) +. ((x_rel *. x_factor) +. x_abs) /. abs (to_real y_f)
|
|
by
|
|
(-. x_rel *. x_factor -. x_abs) /. to_real y_f
|
|
<=. (x_rel *. x_factor +. x_abs) /. abs (to_real y_f)
|
|
};
|
|
div_order_compat2 (x -. x_rel *. x_factor -. x_abs) (to_real x_f) (to_real y_f);
|
|
in ()
|
|
end
|
|
|
|
(** {5 Double propagation lemmas} *)
|
|
module UDoubleLemmas
|
|
use real.RealInfix
|
|
use real.FromInt
|
|
use real.Abs
|
|
use UDouble
|
|
|
|
let lemma uadd_double_error_propagation (x_f y_f r : udouble) (x x_factor x_rel x_abs y y_factor y_rel y_abs : real)
|
|
requires {
|
|
abs (to_real x_f -. x) <=. x_rel *. x_factor +. x_abs
|
|
}
|
|
requires {
|
|
abs (to_real y_f -. y) <=. y_rel *. y_factor +. y_abs
|
|
}
|
|
requires { abs x <=. x_factor }
|
|
requires { abs y <=. y_factor }
|
|
(* TODO: Use (0 <=. x_rel \/ (x_factor = 0 /\ x_abs = 0)), same for y. *)
|
|
requires { 0. <=. x_rel }
|
|
requires { 0. <=. y_rel }
|
|
requires { 0. <=. x_abs }
|
|
requires { 0. <=. y_abs }
|
|
requires { r = x_f ++. y_f }
|
|
ensures {
|
|
abs (to_real r -. (x +. y)) <=.
|
|
(x_rel +. y_rel +. eps) *. (x_factor +. y_factor)
|
|
+. ((1. +. eps +. y_rel) *. x_abs +. (1. +. eps +. x_rel) *. y_abs)
|
|
}
|
|
=
|
|
let ghost delta = abs (to_real (x_f ++. y_f) -. (to_real x_f +. to_real y_f)) in
|
|
assert {
|
|
0. <=. x_rel /\ 0. <=. y_rel ->
|
|
delta <=.
|
|
(eps +. y_rel) *. x_factor +. (eps +. x_rel) *. y_factor
|
|
+. (x_rel +. eps) *. y_abs +. (y_rel +. eps) *. x_abs
|
|
by
|
|
(delta <=. x_factor *. x_rel +. x_abs +. x_factor
|
|
so
|
|
x_factor +. x_abs <=. eps *. (y_factor +. y_abs) ->
|
|
delta <=. (eps +. x_rel) *. y_factor
|
|
+. (eps +. y_rel) *. x_factor
|
|
+. (y_rel +. eps) *. x_abs +. (x_rel +. eps) *. y_abs
|
|
by
|
|
delta <=. eps *. (y_factor +. y_abs) *. x_rel
|
|
+. (eps *. (y_factor +. y_abs)))
|
|
/\
|
|
(delta <=. y_factor *. y_rel +. y_abs +. y_factor
|
|
so
|
|
abs y_factor +. y_abs <=. eps *. (x_factor +. x_abs) ->
|
|
delta <=. (eps +. y_rel) *. x_factor
|
|
+. (eps +. x_rel) *. y_factor
|
|
+. (x_rel +. eps) *. y_abs +. (y_rel +. eps) *. x_abs
|
|
by
|
|
delta <=. eps *. (x_factor +. x_abs) *. y_rel
|
|
+. (eps *. (x_factor +. x_abs)))
|
|
/\
|
|
(
|
|
(eps *. (x_factor +. x_abs) <. abs y_factor +. y_abs /\
|
|
eps *. (y_factor +. y_abs) <. abs x_factor +. x_abs) ->
|
|
(delta <=.
|
|
(eps +. y_rel) *. x_factor +. (eps +. x_rel) *. y_factor
|
|
+. (x_rel +. eps) *. y_abs +. (y_rel +. eps) *. x_abs
|
|
by
|
|
abs (to_real x_f +. to_real y_f) <=.
|
|
abs (to_real x_f -. x) +. x_factor +. (abs (to_real y_f -. y) +. y_factor)
|
|
so
|
|
x_factor *. x_rel <=. (y_factor +. y_abs) /. eps *. x_rel /\
|
|
y_factor *. y_rel <=. (x_factor +. x_abs) /. eps *. y_rel))
|
|
}
|
|
|
|
let lemma usub_double_error_propagation (x_f y_f r : udouble) (x x_factor x_rel x_abs y y_factor y_rel y_abs : real)
|
|
requires {
|
|
abs (to_real x_f -. x) <=. x_rel *. x_factor +. x_abs
|
|
}
|
|
requires {
|
|
abs (to_real y_f -. y) <=. y_rel *. y_factor +. y_abs
|
|
}
|
|
requires { abs x <=. x_factor }
|
|
requires { abs y <=. y_factor }
|
|
requires { 0. <=. x_abs }
|
|
requires { 0. <=. y_abs }
|
|
requires { 0. <=. x_rel }
|
|
requires { 0. <=. y_rel }
|
|
requires { r = x_f --. y_f }
|
|
ensures {
|
|
abs (to_real r -. (x -. y))
|
|
<=. (x_rel +. y_rel +. eps) *. (x_factor +. y_factor)
|
|
+. ((1. +. eps +. y_rel) *. x_abs +. (1. +. eps +. x_rel) *. y_abs)
|
|
}
|
|
= uadd_double_error_propagation x_f (--. y_f) r x x_factor x_rel x_abs (-. y) y_factor y_rel y_abs
|
|
|
|
use HelperLemmas
|
|
|
|
let lemma umul_double_error_propagation (x_f y_f r : udouble) (x x_factor x_rel x_abs y y_factor y_rel y_abs : real)
|
|
requires {
|
|
abs (to_real x_f -. x) <=. x_rel *. x_factor +. x_abs
|
|
}
|
|
requires {
|
|
abs (to_real y_f -. y) <=. y_rel *. y_factor +. y_abs
|
|
}
|
|
requires { abs x <=. x_factor }
|
|
requires { abs y <=. y_factor }
|
|
requires { 0. <=. x_rel }
|
|
requires { 0. <=. y_rel }
|
|
requires { 0. <=. x_abs }
|
|
requires { 0. <=. y_abs }
|
|
requires { r = x_f **. y_f }
|
|
ensures {
|
|
abs (to_real r -. (x *. y)) <=.
|
|
(eps +. (x_rel +. y_rel +. x_rel *. y_rel) *. (1. +. eps)) *. (x_factor *. y_factor)
|
|
+. (((y_abs +. y_abs *. x_rel) *. x_factor
|
|
+. (x_abs +. x_abs *. y_rel) *. y_factor
|
|
+. x_abs *. y_abs) *. (1. +. eps) +. eta)
|
|
}
|
|
=
|
|
assert {
|
|
to_real x_f *. to_real y_f -. abs (to_real x_f *. to_real y_f) *. eps -. eta
|
|
<=. to_real (x_f **. y_f)
|
|
<=. to_real x_f *. to_real y_f +. abs (to_real x_f *. to_real y_f) *. eps +. eta
|
|
};
|
|
assert { abs (x *. y) <=. x_factor *. y_factor by
|
|
abs x *. abs y <=. x_factor *. abs y = abs y *. x_factor <=. y_factor *. x_factor };
|
|
mult_err_combine (to_real x_f) x x_factor x_rel x_abs (to_real y_f) y y_factor y_rel y_abs
|
|
|
|
use real.ExpLog
|
|
|
|
let lemma log_double_error_propagation (logx_f x_f : udouble)
|
|
(x_exact x_factor log_rel log_abs x_rel x_abs : real)
|
|
requires { abs (to_real x_f -. x_exact) <=. x_rel *. x_factor +. x_abs }
|
|
requires {
|
|
abs (to_real logx_f -. log(to_real x_f))
|
|
<=. log_rel *. abs (log (to_real x_f)) +. log_abs
|
|
}
|
|
requires { 0. <. x_exact <=. x_factor }
|
|
requires { 0. <. (x_exact -. x_rel *. x_factor -. x_abs) }
|
|
requires { 0. <=. log_rel }
|
|
ensures {
|
|
abs (to_real logx_f -. log (x_exact))
|
|
<=. log_rel *. abs (log (x_exact)) +.
|
|
(-. log (1. -. ((x_rel *. x_factor +. x_abs) /. x_exact)) *. (1. +. log_rel)
|
|
+. log_abs)
|
|
}
|
|
=
|
|
log_approx_err x_exact (to_real x_f) x_factor x_rel x_abs;
|
|
assert {
|
|
abs (log (to_real x_f)) *. log_rel
|
|
<=. (abs (log (x_exact)) -. log (1.0 -. (((x_rel *. x_factor) +. x_abs) /. x_exact))) *. log_rel
|
|
}
|
|
|
|
let lemma log2_double_error_propagation (log2x_f x_f : udouble)
|
|
(x_exact x_factor log_rel log_abs x_rel x_abs : real)
|
|
requires { abs (to_real x_f -. x_exact) <=. x_rel *. x_factor +. x_abs }
|
|
requires {
|
|
abs (to_real log2x_f -. log2(to_real x_f))
|
|
<=. log_rel *. abs (log2 (to_real x_f)) +. log_abs
|
|
}
|
|
requires { 0. <. x_exact <=. x_factor }
|
|
requires { 0. <. (x_exact -. x_rel *. x_factor -. x_abs) }
|
|
requires { 0. <=. log_rel }
|
|
ensures {
|
|
abs (to_real log2x_f -. log2 (x_exact))
|
|
<=. log_rel *. abs (log2 (x_exact)) +.
|
|
(-. log2 (1. -. ((x_rel *. x_factor +. x_abs) /. x_exact)) *. (1. +. log_rel)
|
|
+. log_abs)
|
|
}
|
|
=
|
|
log2_approx_err x_exact (to_real x_f) x_factor x_rel x_abs;
|
|
assert {
|
|
abs (log2 (to_real x_f)) *. log_rel
|
|
<=. (abs (log2 (x_exact)) -. log2 (1.0 -. (((x_rel *. x_factor) +. x_abs) /. x_exact))) *. log_rel
|
|
}
|
|
|
|
let lemma log10_double_error_propagation (log10x_f x_f : udouble)
|
|
(x_exact x_factor log_rel log_abs x_rel x_abs : real)
|
|
requires { abs (to_real x_f -. x_exact) <=. x_rel *. x_factor +. x_abs }
|
|
requires {
|
|
abs (to_real log10x_f -. log10(to_real x_f))
|
|
<=. log_rel *. abs (log10 (to_real x_f)) +. log_abs
|
|
}
|
|
requires { 0. <. x_exact <=. x_factor }
|
|
requires { 0. <. (x_exact -. x_rel *. x_factor -. x_abs) }
|
|
requires { 0. <=. log_rel }
|
|
ensures {
|
|
abs (to_real log10x_f -. log10 (x_exact))
|
|
<=. log_rel *. abs (log10 (x_exact)) +.
|
|
(-. log10 (1. -. ((x_rel *. x_factor +. x_abs) /. x_exact)) *. (1. +. log_rel)
|
|
+. log_abs)
|
|
}
|
|
=
|
|
log10_approx_err x_exact (to_real x_f) x_factor x_rel x_abs;
|
|
assert {
|
|
abs (log10 (to_real x_f)) *. log_rel
|
|
<=. (abs (log10 (x_exact)) -. log10 (1.0 -. (((x_rel *. x_factor) +. x_abs) /. x_exact))) *. log_rel
|
|
}
|
|
|
|
let lemma exp_double_error_propagation (expx_f x_f : udouble)
|
|
(x_exact x_factor exp_rel exp_abs x_rel x_abs : real)
|
|
requires { abs (to_real x_f -. x_exact) <=. x_rel *. x_factor +. x_abs }
|
|
requires {
|
|
abs (to_real expx_f -. exp(to_real x_f))
|
|
<=. exp_rel *. exp (to_real x_f) +. exp_abs
|
|
}
|
|
requires { x_exact <=. x_factor }
|
|
requires { 0. <=. exp_rel <=. 1. }
|
|
ensures {
|
|
abs (to_real expx_f -. exp (x_exact))
|
|
<=. (exp_rel +. (exp(x_rel *. x_factor +. x_abs) -. 1.) *. (1. +. exp_rel)) *. exp(x_exact)
|
|
+. exp_abs
|
|
}
|
|
=
|
|
exp_approx_err x_exact (to_real x_f) x_factor x_rel x_abs;
|
|
assert {
|
|
exp x_exact *. (1. -. exp_rel) -.
|
|
exp x_exact *. (exp (x_rel *. x_factor +. x_abs) -. 1.) *. (1. -. exp_rel)
|
|
-. exp_abs
|
|
<=. to_real expx_f
|
|
by
|
|
(exp x_exact -. exp x_exact *. (exp (x_rel *. x_factor +. x_abs) -. 1.))
|
|
*. (1. -. exp_rel) -. exp_abs
|
|
<=. exp (to_real x_f) *. (1. -. exp_rel) -. exp_abs
|
|
<=. to_real expx_f
|
|
};
|
|
assert {
|
|
to_real expx_f <=. (exp(x_exact) +. exp(x_exact)*.(exp(x_rel *. x_factor +. x_abs) -. 1.))*. (1. +. exp_rel) +. exp_abs
|
|
by
|
|
to_real expx_f <=. exp(to_real x_f) *. (1. +. exp_rel) +. exp_abs
|
|
};
|
|
|
|
|
|
use real.Trigonometry
|
|
|
|
let lemma sin_double_error_propagation (sinx_f x_f : udouble)
|
|
(x_exact x_factor sin_rel sin_abs x_rel x_abs : real)
|
|
requires { abs (to_real x_f -. x_exact) <=. x_rel *. x_factor +. x_abs }
|
|
requires {
|
|
abs (to_real sinx_f -. sin(to_real x_f))
|
|
<=. sin_rel *. abs (sin (to_real x_f)) +. sin_abs
|
|
}
|
|
requires { x_exact <=. x_factor }
|
|
requires { 0. <=. sin_rel }
|
|
ensures {
|
|
abs (to_real sinx_f -. sin (x_exact))
|
|
<=. sin_rel *. abs(sin(x_exact))
|
|
+. (((x_rel *. x_factor +. x_abs) *. (1. +. sin_rel)) +. sin_abs)
|
|
}
|
|
=
|
|
assert {
|
|
abs (sin (to_real x_f)) *. sin_rel
|
|
<=. (abs (sin x_exact) +. (x_rel *. x_factor +. x_abs)) *. sin_rel
|
|
}
|
|
|
|
let lemma cos_double_error_propagation (cosx_f x_f : udouble)
|
|
(x_exact x_factor cos_rel cos_abs x_rel x_abs : real)
|
|
requires { abs (to_real x_f -. x_exact) <=. x_rel *. x_factor +. x_abs }
|
|
requires {
|
|
abs (to_real cosx_f -. cos(to_real x_f))
|
|
<=. cos_rel *. abs (cos (to_real x_f)) +. cos_abs
|
|
}
|
|
requires { x_exact <=. x_factor }
|
|
requires { 0. <=. cos_rel }
|
|
ensures {
|
|
abs (to_real cosx_f -. cos (x_exact))
|
|
<=. cos_rel *. abs(cos(x_exact))
|
|
+. (((x_rel *. x_factor +. x_abs) *. (1. +. cos_rel)) +. cos_abs)
|
|
}
|
|
=
|
|
assert {
|
|
abs (cos (to_real x_f)) *. cos_rel
|
|
<=. (abs (cos x_exact) +. (x_rel *. x_factor +. x_abs)) *. cos_rel
|
|
}
|
|
|
|
use real.Sum
|
|
use int.Int
|
|
use real.FromInt
|
|
|
|
function real_fun (f:int -> udouble) : int -> real = fun i -> to_real (f i)
|
|
|
|
let lemma sum_double_error_propagation (x : udouble)
|
|
(f : int -> udouble) (f_exact f_factor f_factor' : int -> real) (n:int)
|
|
(sum_rel sum_abs f_rel f_abs : real)
|
|
requires {
|
|
forall i. 0 <= i < n ->
|
|
abs ((real_fun f) i -. f_exact i) <=. f_rel *. f_factor i +. f_abs
|
|
}
|
|
requires {
|
|
forall i. 0 <= i < n ->
|
|
f_factor i -. f_rel *. f_factor i -. f_abs <=. f_factor' i <=. f_factor i +. f_rel *. f_factor i +. f_abs
|
|
}
|
|
requires {
|
|
abs (to_real x -. (sum (real_fun f) 0 n))
|
|
<=. sum_rel *. (sum f_factor' 0 n) +. sum_abs
|
|
}
|
|
requires { 0. <=. sum_rel }
|
|
requires { 0 <= n }
|
|
ensures {
|
|
abs (to_real x -. sum f_exact 0 n)
|
|
<=. (f_rel +. (sum_rel *. (1. +. f_rel))) *. sum f_factor 0 n +.
|
|
((f_abs *. from_int n *.(1. +. sum_rel)) +. sum_abs)
|
|
}
|
|
=
|
|
sum_approx_err f_rel f_abs (real_fun f) f_exact f_factor 0 n;
|
|
sum_approx_err f_rel f_abs f_factor' f_factor f_factor 0 n;
|
|
assert {
|
|
sum_rel *. sum f_factor' 0 n <=.
|
|
sum_rel *. (sum f_factor 0 n +. ((f_rel *. sum f_factor 0 n) +. (f_abs *. from_int n)))
|
|
}
|
|
|
|
(* We don't have an error on y_f because in practice we won't have an exact division with an approximate divisor *)
|
|
let lemma udiv_exact_single_error_propagation (x_f y_f r: udouble) (x x_factor x_rel x_abs : real)
|
|
requires {
|
|
abs (to_real x_f -. x) <=. x_rel *. x_factor +. x_abs
|
|
}
|
|
requires { abs x <=. x_factor }
|
|
requires { 0. <=. x_rel }
|
|
requires { 0. <=. x_abs }
|
|
requires { 0. <> to_real y_f }
|
|
requires { is_exact udiv x_f y_f }
|
|
requires { r = x_f ///. y_f }
|
|
ensures {
|
|
abs (to_real r -. (x /. (to_real y_f))) <=.
|
|
x_rel *. (x_factor /. abs (to_real y_f)) +. ((x_abs /. abs (to_real y_f)) +. eta)
|
|
}
|
|
=
|
|
let lemma y_f_pos ()
|
|
requires { 0. <. to_real y_f }
|
|
ensures {
|
|
abs (to_real r -. (x /. (to_real y_f))) <=.
|
|
x_rel *. (x_factor /. to_real y_f) +. ((x_abs /. to_real y_f) +. eta)
|
|
}
|
|
=
|
|
div_order_compat (to_real x_f) (x +. x_rel *. x_factor +. x_abs) (to_real y_f);
|
|
div_order_compat (x -. x_rel *. x_factor -. x_abs) (to_real x_f) (to_real y_f) in
|
|
let lemma y_f_neg ()
|
|
requires { to_real y_f <. 0. }
|
|
ensures {
|
|
abs (to_real r -. (x /. (to_real y_f))) <=.
|
|
x_rel *. (x_factor /. abs (to_real y_f)) +. ((x_abs /. abs (to_real y_f)) +. eta)
|
|
}
|
|
=
|
|
div_order_compat2 (to_real x_f) (x +. x_rel *. x_factor +. x_abs) (to_real y_f);
|
|
(* TODO: Prove this somehow *)
|
|
assert {
|
|
forall x y. x /. y <=. abs x /. abs y
|
|
};
|
|
assert {
|
|
(x -. x_rel *. x_factor -. x_abs) /. to_real y_f
|
|
<=. x /. (to_real y_f) +. ((x_rel *. x_factor) +. x_abs) /. abs (to_real y_f)
|
|
by
|
|
(-. x_rel *. x_factor -. x_abs) /. to_real y_f
|
|
<=. (x_rel *. x_factor +. x_abs) /. abs (to_real y_f)
|
|
};
|
|
div_order_compat2 (x -. x_rel *. x_factor -. x_abs) (to_real x_f) (to_real y_f);
|
|
in ()
|
|
end
|