From e26fc7430b2e56e04f434ecb992e1c49b4fcd3da Mon Sep 17 00:00:00 2001 From: Michael Tautschnig Date: Tue, 6 Feb 2024 13:29:19 +0000 Subject: [PATCH] C library: model __builtin_powi{,f,l} These GCC built-ins are modelled by implementing pow{,f,l} specialised to integral exponents (where multiple failure cases cannot occur). --- .../cbmc-library/__builtin_powi-01/main.c | 11 + .../cbmc-library/__builtin_powi-01/test.desc | 8 + .../cbmc-library/__builtin_powif-01/main.c | 11 + .../cbmc-library/__builtin_powif-01/test.desc | 8 + .../cbmc-library/__builtin_powil-01/main.c | 11 + .../cbmc-library/__builtin_powil-01/test.desc | 8 + .../loop_contracts_synthesis_04/test.desc | 2 +- src/ansi-c/library/math.c | 361 ++++++++++++++++++ 8 files changed, 419 insertions(+), 1 deletion(-) create mode 100644 regression/cbmc-library/__builtin_powi-01/main.c create mode 100644 regression/cbmc-library/__builtin_powi-01/test.desc create mode 100644 regression/cbmc-library/__builtin_powif-01/main.c create mode 100644 regression/cbmc-library/__builtin_powif-01/test.desc create mode 100644 regression/cbmc-library/__builtin_powil-01/main.c create mode 100644 regression/cbmc-library/__builtin_powil-01/test.desc diff --git a/regression/cbmc-library/__builtin_powi-01/main.c b/regression/cbmc-library/__builtin_powi-01/main.c new file mode 100644 index 00000000000..93b086ce1a6 --- /dev/null +++ b/regression/cbmc-library/__builtin_powi-01/main.c @@ -0,0 +1,11 @@ +#include +#include + +double __builtin_powi(double, int); + +int main() +{ + double four = __builtin_powi(2.0, 2); + assert(four > 3.999 && four < 4.345); + return 0; +} diff --git a/regression/cbmc-library/__builtin_powi-01/test.desc b/regression/cbmc-library/__builtin_powi-01/test.desc new file mode 100644 index 00000000000..3510d48c5c6 --- /dev/null +++ b/regression/cbmc-library/__builtin_powi-01/test.desc @@ -0,0 +1,8 @@ +CORE +main.c +--float-overflow-check --nan-check +^EXIT=0$ +^SIGNAL=0$ +^VERIFICATION SUCCESSFUL$ +-- +^warning: ignoring diff --git a/regression/cbmc-library/__builtin_powif-01/main.c b/regression/cbmc-library/__builtin_powif-01/main.c new file mode 100644 index 00000000000..9eaaba52adb --- /dev/null +++ b/regression/cbmc-library/__builtin_powif-01/main.c @@ -0,0 +1,11 @@ +#include +#include + +float __builtin_powif(float, int); + +int main() +{ + float four = __builtin_powif(2.0f, 2); + assert(four > 3.999f && four < 4.345f); + return 0; +} diff --git a/regression/cbmc-library/__builtin_powif-01/test.desc b/regression/cbmc-library/__builtin_powif-01/test.desc new file mode 100644 index 00000000000..3510d48c5c6 --- /dev/null +++ b/regression/cbmc-library/__builtin_powif-01/test.desc @@ -0,0 +1,8 @@ +CORE +main.c +--float-overflow-check --nan-check +^EXIT=0$ +^SIGNAL=0$ +^VERIFICATION SUCCESSFUL$ +-- +^warning: ignoring diff --git a/regression/cbmc-library/__builtin_powil-01/main.c b/regression/cbmc-library/__builtin_powil-01/main.c new file mode 100644 index 00000000000..cea34ac74e1 --- /dev/null +++ b/regression/cbmc-library/__builtin_powil-01/main.c @@ -0,0 +1,11 @@ +#include +#include + +long double __builtin_powil(long double, int); + +int main() +{ + long double four = __builtin_powil(2.0l, 2); + assert(four > 3.999l && four < 4.345l); + return 0; +} diff --git a/regression/cbmc-library/__builtin_powil-01/test.desc b/regression/cbmc-library/__builtin_powil-01/test.desc new file mode 100644 index 00000000000..3510d48c5c6 --- /dev/null +++ b/regression/cbmc-library/__builtin_powil-01/test.desc @@ -0,0 +1,8 @@ +CORE +main.c +--float-overflow-check --nan-check +^EXIT=0$ +^SIGNAL=0$ +^VERIFICATION SUCCESSFUL$ +-- +^warning: ignoring diff --git a/regression/goto-synthesizer/loop_contracts_synthesis_04/test.desc b/regression/goto-synthesizer/loop_contracts_synthesis_04/test.desc index 9ebb1913e74..2497f2a2345 100644 --- a/regression/goto-synthesizer/loop_contracts_synthesis_04/test.desc +++ b/regression/goto-synthesizer/loop_contracts_synthesis_04/test.desc @@ -3,7 +3,7 @@ main.c --pointer-check _ --no-malloc-may-fail --verbosity 9 ^EXIT=0$ ^SIGNAL=0$ -Quick filter\: 6.* out of 67 candidates were filtered out\. +Quick filter\: [1-9]\d* out of 67 candidates were filtered out\. ^\[main.pointer\_dereference.\d+\] .* SUCCESS$ ^VERIFICATION SUCCESSFUL$ -- diff --git a/src/ansi-c/library/math.c b/src/ansi-c/library/math.c index 32a460297e5..efcf9fc4244 100644 --- a/src/ansi-c/library/math.c +++ b/src/ansi-c/library/math.c @@ -3998,3 +3998,364 @@ long double fmal(long double x, long double y, long double z) return x_times_y + z; #endif } + +/* FUNCTION: __builtin_powi */ + +#ifndef __CPROVER_MATH_H_INCLUDED +# include +# define __CPROVER_MATH_H_INCLUDED +#endif + +#ifndef __CPROVER_STDINT_H_INCLUDED +# include +# define __CPROVER_STDINT_H_INCLUDED +#endif + +#ifndef __CPROVER_ERRNO_H_INCLUDED +# include +# define __CPROVER_ERRNO_H_INCLUDED +#endif + +int32_t __VERIFIER_nondet_int32_t(void); + +double __builtin_inf(void); + +double __builtin_powi(double x, int y) +{ + // see man pow (https://linux.die.net/man/3/pow), specialized for y being an + // integer + if(x == 1.0) + return 1.0; + else if(y == 0) + return 1.0; + else if(fpclassify(x) == FP_ZERO && y > 0) + { + if(y % 2 == 1) + return x; + else + return +0.0; + } + else if(isinf(x) && __CPROVER_signd(x)) + { + if(y < 0) + { + if(-y % 2 == 1) + return -0.0; + else + return +0.0; + } + else + { + if(y % 2 == 1) + return -__builtin_inf(); + else + return __builtin_inf(); + } + } + else if(isinf(x) && !__CPROVER_signd(x)) + { + if(y < 0) + return +0.0; + else + return __builtin_inf(); + } + else if(fpclassify(x) == FP_ZERO && y < 0) + { + errno = ERANGE; +#pragma CPROVER check push +#pragma CPROVER check disable "float-overflow" + if(__CPROVER_signd(x) && -y % 2 == 1) + return -HUGE_VAL; + else + return HUGE_VAL; +#pragma CPROVER check pop + } + else if(isnan(x)) +#pragma CPROVER check push +#pragma CPROVER check disable "div-by-zero" + return 0.0 / 0.0; +#pragma CPROVER check pop + +#ifndef _MSC_VER + _Static_assert +#else + static_assert +#endif + (sizeof(double) == 2 * sizeof(int32_t), + "bit width of double is 2x bit width of int32_t"); + // https://martin.ankerl.com/2007/10/04/optimized-pow-approximation-for-java-and-c-c/ + union + { + double d; + int32_t i[2]; + } u = {x}; + int32_t bias = (1 << 20) * ((1 << 10) - 1); + int32_t exp_c = __VERIFIER_nondet_int32_t(); + __CPROVER_assume(exp_c >= -90253 && exp_c <= 1); +#pragma CPROVER check push +#pragma CPROVER check disable "signed-overflow" +#if !defined(__BYTE_ORDER__) || __BYTE_ORDER__ == __ORDER_LITTLE_ENDIAN__ + double mult_result = (double)(y) * (double)(u.i[1] - (bias + exp_c)); +#else + double mult_result = (double)(y) * (double)(u.i[0] - (bias + exp_c)); +#endif +#pragma CPROVER check pop + if(fabs(mult_result) > (double)(1 << 30)) + { + errno = ERANGE; +#pragma CPROVER check push +#pragma CPROVER check disable "float-overflow" + return y > 0 ? HUGE_VAL : 0.0; +#pragma CPROVER check pop + } +#if !defined(__BYTE_ORDER__) || __BYTE_ORDER__ == __ORDER_LITTLE_ENDIAN__ + u.i[1] = (int32_t)mult_result + (bias + exp_c); + u.i[0] = 0; +#else + u.i[0] = (int32_t)mult_result + (bias + exp_c); + u.i[1] = 0; +#endif + return u.d; +} + +/* FUNCTION: __builtin_powif */ + +#ifndef __CPROVER_MATH_H_INCLUDED +# include +# define __CPROVER_MATH_H_INCLUDED +#endif + +#ifndef __CPROVER_STDINT_H_INCLUDED +# include +# define __CPROVER_STDINT_H_INCLUDED +#endif + +#ifndef __CPROVER_ERRNO_H_INCLUDED +# include +# define __CPROVER_ERRNO_H_INCLUDED +#endif + +int32_t __VERIFIER_nondet_int32_t(void); + +float __builtin_inff(void); + +float __builtin_powif(float x, int y) +{ + // see man pow (https://linux.die.net/man/3/pow), specialized for y being an + // integer + if(x == 1.0f) + return 1.0f; + else if(y == 0) + return 1.0f; + else if(fpclassify(x) == FP_ZERO && y > 0) + { + if(y % 2 == 1) + return x; + else + return +0.0f; + } + else if(isinff(x) && __CPROVER_signf(x)) + { + if(y < 0) + { + if(-y % 2 == 1) + return -0.0f; + else + return +0.0f; + } + else + { + if(y % 2 == 1) + return -__builtin_inff(); + else + return __builtin_inff(); + } + } + else if(isinff(x) && !__CPROVER_signf(x)) + { + if(y < 0) + return +0.0f; + else + return __builtin_inff(); + } + else if(fpclassify(x) == FP_ZERO && y < 0) + { + errno = ERANGE; +#pragma CPROVER check push +#pragma CPROVER check disable "float-overflow" + if(__CPROVER_signf(x) && -y % 2 == 1) + { + return -HUGE_VALF; + } + else + return HUGE_VALF; +#pragma CPROVER check pop + } + else if(isnanf(x)) +#pragma CPROVER check push +#pragma CPROVER check disable "div-by-zero" + return 0.0f / 0.0f; +#pragma CPROVER check pop + +#ifndef _MSC_VER + _Static_assert +#else + static_assert +#endif + (sizeof(float) == sizeof(int32_t), + "bit width of float and int32_t should match"); + union + { + float f; + int32_t i; + } u = {x}; + int32_t bias = (1 << 23) * ((1 << 7) - 1); + int32_t exp_c = __VERIFIER_nondet_int32_t(); + __CPROVER_assume(exp_c >= -722019 && exp_c <= 1); +#pragma CPROVER check push +#pragma CPROVER check disable "signed-overflow" + float mult_result = (float)(y) * (float)(u.i - (bias + exp_c)); +#pragma CPROVER check pop + if(fabsf(mult_result) > (float)(1 << 30)) + { + errno = ERANGE; +#pragma CPROVER check push +#pragma CPROVER check disable "float-overflow" + return y > 0 ? HUGE_VALF : 0.0f; +#pragma CPROVER check pop + } + u.i = (int32_t)mult_result + (bias + exp_c); + return u.f; +} + +/* FUNCTION: __builtin_powil */ + +#ifndef __CPROVER_MATH_H_INCLUDED +# include +# define __CPROVER_MATH_H_INCLUDED +#endif + +#ifndef __CPROVER_FLOAT_H_INCLUDED +# include +# define __CPROVER_FLOAT_H_INCLUDED +#endif + +#ifndef __CPROVER_STDINT_H_INCLUDED +# include +# define __CPROVER_STDINT_H_INCLUDED +#endif + +#ifndef __CPROVER_ERRNO_H_INCLUDED +# include +# define __CPROVER_ERRNO_H_INCLUDED +#endif + +int32_t __VERIFIER_nondet_int32_t(void); + +long double __builtin_infl(void); +double __builtin_powi(double, int); + +long double __builtin_powil(long double x, int y) +{ + // see man pow (https://linux.die.net/man/3/pow), specialized for y being an + // integer + if(x == 1.0l) + return 1.0l; + else if(y == 0) + return 1.0l; + else if(fpclassify(x) == FP_ZERO && y > 0) + { + if(y % 2 == 1) + return x; + else + return +0.0l; + } + else if(isinf(x) && __CPROVER_signld(x)) + { + if(y < 0) + { + if(-y % 2 == 1) + return -0.0l; + else + return +0.0l; + } + else + { + if(y % 2 == 1) + return -__builtin_infl(); + else + return __builtin_infl(); + } + } + else if(isinf(x) && !__CPROVER_signld(x)) + { + if(y < 0) + return +0.0f; + else + return __builtin_infl(); + } + else if(fpclassify(x) == FP_ZERO && y < 0) + { + errno = ERANGE; +#pragma CPROVER check push +#pragma CPROVER check disable "float-overflow" + if(__CPROVER_signld(x) && -y % 2 == 1) + { + return -HUGE_VALL; + } + else + return HUGE_VALL; +#pragma CPROVER check pop + } + else if(isnan(x)) +#pragma CPROVER check push +#pragma CPROVER check disable "div-by-zero" + return 0.0l / 0.0l; +#pragma CPROVER check pop + +#if LDBL_MAX_EXP == DBL_MAX_EXP + return __builtin_powi(x, y); +#else +# ifndef _MSC_VER + _Static_assert +# else + static_assert +# endif + (sizeof(long double) % sizeof(int32_t) == 0, + "bit width of long double is a multiple of bit width of int32_t"); + union U + { + long double l; + int32_t i[sizeof(long double) / sizeof(int32_t)]; + } u = {x}; +# if !defined(__BYTE_ORDER__) || __BYTE_ORDER__ == __ORDER_LITTLE_ENDIAN__ + int32_t exponent = u.i[sizeof(long double) / sizeof(int32_t) - 1]; +# else + int32_t exponent = u.i[0]; +# endif + int32_t bias = (1 << 16) * ((1 << 14) - 1); + int32_t exp_c = __VERIFIER_nondet_int32_t(); + __CPROVER_assume(exp_c >= -5641 && exp_c <= 1); +# pragma CPROVER check push +# pragma CPROVER check disable "signed-overflow" + long double mult_result = + (long double)y * (long double)(exponent - (bias + exp_c)); +# pragma CPROVER check pop + if(fabsl(mult_result) > (long double)(1 << 30)) + { + errno = ERANGE; +# pragma CPROVER check push +# pragma CPROVER check disable "float-overflow" + return y > 0 ? HUGE_VALL : 0.0l; +# pragma CPROVER check pop + } + int32_t result = (int32_t)mult_result + (bias + exp_c); + union U result_u = {.i = {0}}; +# if !defined(__BYTE_ORDER__) || __BYTE_ORDER__ == __ORDER_LITTLE_ENDIAN__ + result_u.i[sizeof(long double) / sizeof(int32_t) - 1] = result; +# else + result_u.i[0] = result; +# endif + return result_u.l; +#endif +}