diff --git a/gobview b/gobview index b5ab4dbd7b..74e2958708 160000 --- a/gobview +++ b/gobview @@ -1 +1 @@ -Subproject commit b5ab4dbd7b1b1dc4f0bf2a912012cb1879f54903 +Subproject commit 74e295870848c66885fa446970b0c59617a242a7 diff --git a/src/cdomain/value/cdomains/floatDomain.ml b/src/cdomain/value/cdomains/floatDomain.ml index c22b6dfa4d..503feba88e 100644 --- a/src/cdomain/value/cdomains/floatDomain.ml +++ b/src/cdomain/value/cdomains/floatDomain.ml @@ -5,6 +5,23 @@ open FloatOps exception ArithmeticOnFloatBot of string +(** Define records that hold mutable variables representing different Configuration values. + * These values are used to keep track of whether or not the corresponding Config values are en-/disabled *) +type ana_float_config_values = { + mutable evaluate_math_functions : bool option; +} + +let ana_float_config: ana_float_config_values = { + evaluate_math_functions = None; +} + +let get_evaluate_math_functions () = + if ana_float_config.evaluate_math_functions = None then + ana_float_config.evaluate_math_functions <- Some (GobConfig.get_bool "ana.float.evaluate_math_functions"); + Option.get ana_float_config.evaluate_math_functions + +let reset_lazy () = ana_float_config.evaluate_math_functions <- None + module type FloatArith = sig type t @@ -273,12 +290,12 @@ module FloatIntervalImpl(Float_t : CFloatType) = struct | _ -> Bot (** [widen x y] assumes [leq x y]. Solvers guarantee this by calling [widen old (join old new)]. *) - let widen v1 v2 = (**TODO: support 'threshold_widening' option *) + let widen v1 v2 = (* TODO: support 'threshold_widening' option *) match v1, v2 with | Top, _ | _, Top -> Top | Bot, v | v, Bot -> v | Interval (l1, h1), Interval (l2, h2) -> - (**If we widen and we know that neither interval contains +-inf or nan, it is ok to widen only to +-max_float, + (* If we widen and we know that neither interval contains +-inf or nan, it is ok to widen only to +-max_float, because a widening with +-inf/nan will always result in the case above -> Top *) let low = if l1 <= l2 then l1 else Float_t.lower_bound in let high = if h1 >= h2 then h1 else Float_t.upper_bound in @@ -289,7 +306,7 @@ module FloatIntervalImpl(Float_t : CFloatType) = struct | _ -> Top let narrow v1 v2 = - match v1, v2 with (**we cannot distinguish between the lower bound beeing -inf or the upper bound beeing inf. Also there is nan *) + match v1, v2 with (* we cannot distinguish between the lower bound beeing -inf or the upper bound beeing inf. Also there is nan *) | Bot, _ | _, Bot -> Bot | Top, _ -> v2 | Interval (l1, h1), Interval (l2, h2) -> @@ -623,15 +640,79 @@ module FloatIntervalImpl(Float_t : CFloatType) = struct else unknown_IInt () - (**it seems strange not to return an explicit 1 for negative numbers, but in c99 signbit is defined as: *) - (**<> *) + (**it seems strange not to return an explicit 1 for negative numbers, but in c99 signbit is defined as: + **<> *) let eval_signbit = function | (_, h) when h < Float_t.zero -> true_nonZero_IInt () | (l, _) when l > Float_t.zero -> false_zero_IInt () - | _ -> unknown_IInt () (**any interval containing zero has to fall in this case, because we do not distinguish between 0. and -0. *) - - (**This Constant overapproximates pi to use as bounds for the return values of trigonometric functions *) - let overapprox_pi = 3.1416 + | _ -> unknown_IInt () (* any interval containing zero has to fall in this case, because we do not distinguish between 0. and -0. *) + + (** returns the max of the given mfun once computed with rounding mode Up and once with rounding mode down*) + let safe_mathfun_up mfun f = max (mfun Down f) (mfun Up f) + (** returns the min of the given mfun once computed with rounding mode Up and once with rounding mode down*) + let safe_mathfun_down mfun f = min (mfun Down f) (mfun Up f) + + (** This function does two things: + ** 1. projects l and h onto the interval [0, k*pi] (for k = 2 this is the phase length of sin/cos, for k = 1 it is the phase length of tan) + ** 2. compresses/transforms the interval [0, k*pi] to the interval [0, 1] to ease further computations + ** i.e. the function computes dist = distance, l'' = (l/(k*pi)) - floor(l/(k*pi)), h'' = (h/(k*pi)) - floor(h/(k*pi))*) + let project_and_compress l h k = + let ft_over_kpi = (Float_t.mul Up (Float_t.of_float Up k) Float_t.pi) in + let ft_under_kpi = (Float_t.mul Down (Float_t.of_float Down k) Float_t.pi) in + let l' = + if l >= Float_t.zero then (Float_t.div Down l ft_over_kpi) + else (Float_t.div Down l ft_under_kpi) + in + let h' = + if h >= Float_t.zero then (Float_t.div Up h ft_under_kpi) + else (Float_t.div Up h ft_over_kpi) + in + let dist = (Float_t.sub Up h' l') in + let l'' = + if l' >= Float_t.zero then + Float_t.sub Down l' (Float_t.of_float Up (Z.to_float (Float_t.to_big_int l'))) + else + Float_t.sub Down l' (Float_t.of_float Up (Z.to_float (Z.pred (Float_t.to_big_int l')))) + in + let h'' = + if h' >= Float_t.zero then + Float_t.sub Up h' (Float_t.of_float Down (Z.to_float (Float_t.to_big_int h'))) + else + Float_t.sub Up h' (Float_t.of_float Down (Z.to_float (Z.pred (Float_t.to_big_int h')))) + in + (dist, l'', h'') + + let eval_cos_cfun l h = + let (dist, l'', h'') = project_and_compress l h 2. in + if Messages.tracing then Messages.trace "CstubsTrig" "cos: dist %s; l'' %s; h'' %s\n" (Float_t.to_string dist) (Float_t.to_string l'') (Float_t.to_string h''); + if (dist <= Float_t.of_float Down 0.5) && (h'' <= Float_t.of_float Down 0.5) && (l'' <= h'') then + (* case: monotonic decreasing interval*) + Interval (safe_mathfun_down Float_t.cos h, safe_mathfun_up Float_t.cos l) + else if (dist <= Float_t.of_float Down 0.5) && (l'' >= Float_t.of_float Up 0.5) && (l'' <= h'') then + (* case: monotonic increasing interval*) + Interval (safe_mathfun_down Float_t.cos l, safe_mathfun_up Float_t.cos h) + else if (dist <= Float_t.of_float Down 1.) && (l'' <= h'') then + (* case: contains at most one minimum*) + Interval (Float_t.of_float Down (-.1.), max (safe_mathfun_up Float_t.cos l) (safe_mathfun_up Float_t.cos h)) + else if (dist <= Float_t.of_float Down 1.) && (l'' >= Float_t.of_float Up 0.5) && (h'' <= Float_t.of_float Down 0.5) then + (* case: contains at most one maximum*) + Interval (min (safe_mathfun_down Float_t.cos l) (safe_mathfun_down Float_t.cos h), Float_t.of_float Up 1.) + else + of_interval (-. 1., 1.) + + let eval_sin_cfun l h = + let lcos = Float_t.sub Down l (Float_t.div Up Float_t.pi (Float_t.of_float Down 2.0)) in + let hcos = Float_t.sub Up h (Float_t.div Down Float_t.pi (Float_t.of_float Up 2.0)) in + eval_cos_cfun lcos hcos + + let eval_tan_cfun l h = + let (dist, l'', h'') = project_and_compress l h 1. in + if Messages.tracing then Messages.trace "CstubsTrig" "tan: dist %s; l'' %s; h'' %s\n" (Float_t.to_string dist) (Float_t.to_string l'') (Float_t.to_string h''); + if (dist <= Float_t.of_float Down 1.) && (Bool.not ((l'' <= Float_t.of_float Up 0.5) && (h'' >= Float_t.of_float Up 0.5))) then + (* case: monotonic increasing interval*) + Interval (safe_mathfun_down Float_t.tan l, safe_mathfun_up Float_t.tan h) + else + top () let eval_fabs = function | (l, h) when l > Float_t.zero -> Interval (l, h) @@ -656,39 +737,51 @@ module FloatIntervalImpl(Float_t : CFloatType) = struct let eval_acos = function | (l, h) when l = h && l = Float_t.of_float Nearest 1. -> of_const 0. (*acos(1) = 0*) - | (l, h) -> - if l < (Float_t.of_float Down (-.1.)) || h > (Float_t.of_float Up 1.) then - Messages.warn ~category:Messages.Category.Float "Domain error might occur: acos argument might be outside of [-1., 1.]"; - of_interval (0., (overapprox_pi)) (**could be more exact *) + | (l, h) when l < (Float_t.of_float Down (-.1.)) || h > (Float_t.of_float Up 1.) -> + Messages.warn ~category:Messages.Category.Float "Domain error might occur: acos argument might be outside of [-1., 1.]"; + Interval (Float_t.of_float Down 0., Float_t.pi) + | (l, h) when get_evaluate_math_functions () -> + norm @@ Interval (safe_mathfun_down Float_t.acos h, safe_mathfun_up Float_t.acos l) (* acos is monotonic decreasing in [-1, 1]*) + | _ -> Interval (Float_t.of_float Down 0., Float_t.pi) let eval_asin = function | (l, h) when l = h && l = Float_t.zero -> of_const 0. (*asin(0) = 0*) - | (l, h) -> - if l < (Float_t.of_float Down (-.1.)) || h > (Float_t.of_float Up 1.) then - Messages.warn ~category:Messages.Category.Float "Domain error might occur: asin argument might be outside of [-1., 1.]"; - div (of_interval ((-. overapprox_pi), overapprox_pi)) (of_const 2.) (**could be more exact *) + | (l, h) when l < (Float_t.of_float Down (-.1.)) || h > (Float_t.of_float Up 1.) -> + Messages.warn ~category:Messages.Category.Float "Domain error might occur: asin argument might be outside of [-1., 1.]"; + div (Interval (Float_t.neg Float_t.pi, Float_t.pi)) (of_const 2.) + | (l, h) when get_evaluate_math_functions () -> + norm @@ Interval (safe_mathfun_down Float_t.asin l, safe_mathfun_up Float_t.asin h) (* asin is monotonic increasing in [-1, 1]*) + | _ -> div (Interval (Float_t.neg Float_t.pi, Float_t.pi)) (of_const 2.) let eval_atan = function | (l, h) when l = h && l = Float_t.zero -> of_const 0. (*atan(0) = 0*) - | _ -> div (of_interval ((-. overapprox_pi), overapprox_pi)) (of_const 2.) (**could be more exact *) + | (l, h) when get_evaluate_math_functions () -> + norm @@ Interval (safe_mathfun_down Float_t.atan l, safe_mathfun_up Float_t.atan h) (* atan is monotonic increasing*) + | _ -> div (Interval (Float_t.neg Float_t.pi, Float_t.pi)) (of_const 2.) let eval_cos = function | (l, h) when l = h && l = Float_t.zero -> of_const 1. (*cos(0) = 1*) - | _ -> of_interval (-. 1., 1.) (**could be exact for intervals where l=h, or even for Interval intervals *) + | (l, h) when get_evaluate_math_functions () -> + norm @@ eval_cos_cfun l h + | _ -> of_interval (-. 1., 1.) let eval_sin = function | (l, h) when l = h && l = Float_t.zero -> of_const 0. (*sin(0) = 0*) - | _ -> of_interval (-. 1., 1.) (**could be exact for intervals where l=h, or even for some intervals *) + | (l, h) when get_evaluate_math_functions () -> + norm @@ eval_sin_cfun l h + | _ -> of_interval (-. 1., 1.) let eval_tan = function | (l, h) when l = h && l = Float_t.zero -> of_const 0. (*tan(0) = 0*) - | _ -> top () (**could be exact for intervals where l=h, or even for some intervals *) + | (l, h) when get_evaluate_math_functions () -> + norm @@ eval_tan_cfun l h + | _ -> top () let eval_sqrt = function | (l, h) when l = Float_t.zero && h = Float_t.zero -> of_const 0. - | (l, h) when l >= Float_t.zero -> - let low = Float_t.sqrt Down l in - let high = Float_t.sqrt Up h in + | (l, h) when l >= Float_t.zero && get_evaluate_math_functions () -> + let low = safe_mathfun_down Float_t.sqrt l in + let high = safe_mathfun_up Float_t.sqrt h in Interval (low, high) | _ -> top () diff --git a/src/cdomain/value/cdomains/floatDomain.mli b/src/cdomain/value/cdomains/floatDomain.mli index d958e1ee81..b672bbdf34 100644 --- a/src/cdomain/value/cdomains/floatDomain.mli +++ b/src/cdomain/value/cdomains/floatDomain.mli @@ -4,6 +4,8 @@ open GoblintCil exception ArithmeticOnFloatBot of string +val reset_lazy: unit -> unit + module type FloatArith = sig type t diff --git a/src/common/cdomains/floatOps/floatOps.ml b/src/common/cdomains/floatOps/floatOps.ml index a791778857..ecfc6436d4 100644 --- a/src/common/cdomains/floatOps/floatOps.ml +++ b/src/common/cdomains/floatOps/floatOps.ml @@ -13,6 +13,7 @@ module type CFloatType = sig val upper_bound: t val lower_bound: t val smallest : t + val pi : t val of_float: round_mode -> float -> t val to_float: t -> float option @@ -37,6 +38,14 @@ module type CFloatType = sig val mul: round_mode -> t -> t -> t val div: round_mode -> t -> t -> t val sqrt: round_mode -> t -> t + + val acos: round_mode -> t -> t + val asin: round_mode -> t -> t + val atan: round_mode -> t -> t + val cos: round_mode -> t -> t + val sin: round_mode -> t -> t + val tan: round_mode -> t -> t + val atof: round_mode -> string -> t end @@ -56,6 +65,7 @@ module CDouble = struct let upper_bound = Float.max_float let lower_bound = -. Float.max_float let smallest = Float.min_float + let pi = Float.pi let of_float _ x = x let to_float x = Some x @@ -78,6 +88,13 @@ module CDouble = struct external div: round_mode -> t -> t -> t = "div_double" external sqrt: round_mode -> t -> t = "sqrt_double" + external acos: round_mode -> t -> t = "acos_double" + external asin: round_mode -> t -> t = "asin_double" + external atan: round_mode -> t -> t = "atan_double" + external cos: round_mode -> t -> t = "cos_double" + external sin: round_mode -> t -> t = "sin_double" + external tan: round_mode -> t -> t = "tan_double" + external atof: round_mode -> string -> t = "atof_double" end @@ -89,10 +106,12 @@ module CFloat = struct external upper': unit -> float = "max_float" external smallest': unit -> float = "smallest_float" + external pi': unit -> float = "pi_float" let upper_bound = upper' () let lower_bound = -. upper_bound let smallest = smallest' () + let pi = pi' () let to_float x = Some x let to_big_int = big_int_of_float @@ -112,6 +131,13 @@ module CFloat = struct external div: round_mode -> t -> t -> t = "div_float" external sqrt: round_mode -> t -> t = "sqrt_float" + external acos: round_mode -> t -> t = "acos_float" + external asin: round_mode -> t -> t = "asin_float" + external atan: round_mode -> t -> t = "atan_float" + external cos: round_mode -> t -> t = "cos_float" + external sin: round_mode -> t -> t = "sin_float" + external tan: round_mode -> t -> t = "tan_float" + external atof: round_mode -> string -> t = "atof_float" let of_float mode x = add mode zero x diff --git a/src/common/cdomains/floatOps/floatOps.mli b/src/common/cdomains/floatOps/floatOps.mli index cf24f75ed5..7cd74f7e7d 100644 --- a/src/common/cdomains/floatOps/floatOps.mli +++ b/src/common/cdomains/floatOps/floatOps.mli @@ -14,6 +14,7 @@ module type CFloatType = sig val upper_bound: t val lower_bound: t val smallest : t + val pi : t val of_float: round_mode -> float -> t val to_float: t -> float option @@ -39,6 +40,12 @@ module type CFloatType = sig val mul: round_mode -> t -> t -> t val div: round_mode -> t -> t -> t val sqrt: round_mode -> t -> t + val acos: round_mode -> t -> t + val asin: round_mode -> t -> t + val atan: round_mode -> t -> t + val cos: round_mode -> t -> t + val sin: round_mode -> t -> t + val tan: round_mode -> t -> t val atof: round_mode -> string -> t end diff --git a/src/common/cdomains/floatOps/stubs.c b/src/common/cdomains/floatOps/stubs.c index b98b26c1da..213eacba22 100644 --- a/src/common/cdomains/floatOps/stubs.c +++ b/src/common/cdomains/floatOps/stubs.c @@ -1,3 +1,4 @@ +#define _GNU_SOURCE // necessary for M_PI to be defined #include #include #include @@ -54,6 +55,18 @@ static void change_round_mode(int mode) UNARY_OP(sqrt, double, sqrt); UNARY_OP(sqrt, float, sqrtf); +UNARY_OP(acos, double, acos); +UNARY_OP(acos, float, acosf); +UNARY_OP(asin, double, asin); +UNARY_OP(asin, float, asinf); +UNARY_OP(atan, double, atan); +UNARY_OP(atan, float, atanf); +UNARY_OP(cos, double, cos); +UNARY_OP(cos, float, cosf); +UNARY_OP(sin, double, sin); +UNARY_OP(sin, float, sinf); +UNARY_OP(tan, double, tan); +UNARY_OP(tan, float, tanf); #define BINARY_OP(name, type, op) \ CAMLprim value name##_##type(value mode, value x, value y) \ @@ -109,6 +122,8 @@ CAMLprim value atof_float(value mode, value str) // No need to use CAMLreturn because we don't use CAMLparam. } +// These are only given for floats as these operations involve no rounding and their OCaml implementation (Float module) can be used + CAMLprim value max_float(value unit) { // No need to use CAMLparam to keep unit as GC root, @@ -124,3 +139,8 @@ CAMLprim value smallest_float(value unit) return caml_copy_double(FLT_MIN); // No need to use CAMLreturn because we don't use CAMLparam. } + +CAMLprim value pi_float(value unit) +{ + return caml_copy_double(M_PI); +} diff --git a/src/config/options.schema.json b/src/config/options.schema.json index a01a5bdb31..8e4ca51fee 100644 --- a/src/config/options.schema.json +++ b/src/config/options.schema.json @@ -463,6 +463,13 @@ "Use FloatDomain: (float * float) option.", "type": "boolean", "default": false + }, + "evaluate_math_functions": { + "title": "ana.float.evaluate_math_functions", + "description": + "Allow a more precise evaluation of some functions from math.h. Evaluation of functions may differ depending on the implementation of math.h. Caution: For some implementations of functions it is not guaranteed that they behave monotonic where they mathematically should, thus possibly leading to unsoundness.", + "type": "boolean", + "default": false } }, "additionalProperties": false diff --git a/src/util/server.ml b/src/util/server.ml index 31e1036652..6f760e6856 100644 --- a/src/util/server.ml +++ b/src/util/server.ml @@ -284,6 +284,7 @@ let analyze ?(reset=false) (s: t) = InvariantCil.reset_lazy (); WideningThresholds.reset_lazy (); IntDomain.reset_lazy (); + FloatDomain.reset_lazy (); StringDomain.reset_lazy (); PrecisionUtil.reset_lazy (); ApronDomain.reset_lazy (); diff --git a/tests/regression/39-signed-overflows/07-abs-sqrt.c b/tests/regression/39-signed-overflows/07-abs-sqrt.c index 13ed863e51..e1e5462362 100644 --- a/tests/regression/39-signed-overflows/07-abs-sqrt.c +++ b/tests/regression/39-signed-overflows/07-abs-sqrt.c @@ -1,4 +1,4 @@ -//PARAM: --enable ana.int.interval --enable ana.float.interval --set ana.activated[+] tmpSpecial +//PARAM: --enable ana.int.interval --enable ana.float.interval --enable ana.float.evaluate_math_functions --set ana.activated[+] tmpSpecial #include int main() { int data; diff --git a/tests/regression/39-signed-overflows/09-labs-sqrt.c b/tests/regression/39-signed-overflows/09-labs-sqrt.c index 3a4b20a82b..8ec231c879 100644 --- a/tests/regression/39-signed-overflows/09-labs-sqrt.c +++ b/tests/regression/39-signed-overflows/09-labs-sqrt.c @@ -1,4 +1,4 @@ -//PARAM: --enable ana.int.interval --enable ana.float.interval --set ana.activated[+] tmpSpecial +//PARAM: --enable ana.int.interval --enable ana.float.interval --enable ana.float.evaluate_math_functions --set ana.activated[+] tmpSpecial #include int main() { int data; diff --git a/tests/regression/57-floats/23-trig_functions_c-stubs.c b/tests/regression/57-floats/23-trig_functions_c-stubs.c new file mode 100644 index 0000000000..7c764a53e4 --- /dev/null +++ b/tests/regression/57-floats/23-trig_functions_c-stubs.c @@ -0,0 +1,67 @@ +// PARAM: --enable ana.float.interval --enable ana.float.evaluate_math_functions +#include +#include +#include + +int main() +{ + int x; + double s1, s2, s3, c1, c2, sc1, sc2, t1, t2; + //s1: 0.5pi < [2.0, 7.5] < 2.5pi + //s2: 1.5pi < [5.0, 10.5] < 3.5pi + //s3: -0.5pi < [-1.0, -1.0] < 0.5pi + //c1: 0pi < [0.2, 6.1] < 2pi + //c2: pi < [3.3, 9.0] < 3pi + //sc1: 0.5pi < [2.0, 3.0] < pi + //sc2: 4.5pi < [14.5, 15.5] < 5pi + //t1: 0pi-0.1 <= [-0.1, -0.1] <= 0pi+0.1 + //t2: 6pi-0.1 < [18.8, 18.9] < 6pi+0.1 + if (x) { + s1 = 2.0; + s2 = 5.0; + s3 = -1.0; + c1 = 0.2; + c2 = 3.3; + sc1 = 2.0; + sc2 = 14.5; + t1 = -0.1; + t2 = 18.8; + } + else { + s1 = 7.5; + s2 = 10.5; + s3 = 1.0; + c1 = 6.1; + c2 = 9.0; + sc1 = 3.0; + sc2 = 15.5; + t1 = 0.1; + t2 = 18.9; + } + + //acos(x) + __goblint_check(1.4 < acos(0.1) && acos(0.1) < 1.5); // SUCCESS + + //asin(x) + __goblint_check(0.6 < asin(0.6) && asin(0.6) < 0.7); // SUCCESS + + //atan(x) + __goblint_check(0.7 < atan(1.) && atan(1.) < 0.8); // SUCCESS + + //cos(x) + __goblint_check(-1. <= cos(c1) && cos(c1) < 0.99); // SUCCESS + __goblint_check(-0.99 < cos(c2) && cos(c2) <= 1.0); // SUCCESS + __goblint_check(-0.99 < cos(sc1) && cos(sc1) < 0.); // SUCCESS + __goblint_check(-0.99 < cos(sc2) && cos(sc2) < 0.); // SUCCESS + + //sin(x) + __goblint_check(-1. <= sin(s1) && sin(s1) < 0.99); // SUCCESS + __goblint_check(-0.99 < sin(s2) && sin(s2) <= 1.); // SUCCESS + __goblint_check(-0.99 < sin(s3) && sin(s3) <= 0.99); // SUCCESS + __goblint_check(0. < sin(sc1) && sin(sc1) < 0.99); // SUCCESS + __goblint_check(0. < sin(sc2) && sin(sc2) < 0.99); // SUCCESS + + //tan(x) + __goblint_check(-0.11 < tan(t1) && tan(t1) < 0.11 ); // SUCCESS + __goblint_check(-0.1 < tan(t2) && tan(t2) < 0.1 ); // SUCCESS +}