Frama-C news and ideas

A contest and a self-pitying lament

John Regehr was organizing a craziest undefined behavior contest and the results are in. I had an entry in the contest, but I did not win. My entry apparently was too obviously dangerous. As John puts it, “I would have expected a modern C compiler to exploit those undefined behaviors”.

I have explained the problem inherent to constructs of the kind &a - 1 < &a in several different sets of circumstances, to different programmers for different projects, and I was unable to convince any of them that the construct or a variant of it was not harmless and that the code should be changed. But for this contest, the construct was too obviously wrong and too prone to be taken advantage of by an aggressive compiler to win. This is the story of my life, really.

Notable entries

The contest invited C++ entries as well, but fortunately for this blog post, most of the submissions were actually in C.

Runner-up

int main (void) {
  int x = 1;
  while (x) x <<= 1;
  return x;
}

$ bin/toplevel.opt -val t1.c | grep assert
t1.c:3:[kernel] warning: invalid LHS operand for left shift. assert x ≥ 0;

As shown above, Frama-C's value analysis warns for the undefined behavior that made this entry a runner-up in the contest.

Winner #1

enum {N = 32};
int a[N], pfx[N];
void prefix_sum (void)
{
  int i, accum;
  for (i = 0, accum = a[0]; i < N; i++, accum += a[i])
    pfx[i] = accum;
}

$ bin/toplevel.opt -val t2.c -main prefix_sum | grep assert
t2.c:6:[kernel] warning: accessing out of bounds index [1..32]. assert i < N;

In this second example, you need to tell Frama-C to analyze the function prefix_sum. Otherwise, it won't find the undefined behavior: the value analysis is not the sort of heuristic bug-finding tool that detects patterns in code regardless of calling contexts. For precision and soundness, Frama-C's value analysis tries to determine whether prefix_sum is safe in the context(s) in which it is actually called during execution. Without option -main prefix_sum, it looks like it is not called at all, and therefore safe. Indeed, if you try to compile and run the program, the linker complains about a missing main() function, meaning that the code cannot do harm.

Having to use option -main and others to describe your intentions is the sort of quirk you can pick up by going through a tutorial or the documentation.

Winner #2

#include <stdio.h>
#include <stdlib.h>
 
int main() {
  int *p = (int*)malloc(sizeof(int));
  int *q = (int*)realloc(p, sizeof(int));
  *p = 1;
  *q = 2;
  if (p == q)
    printf("%d %d\n", *p, *q);
}

This example involves dynamic allocation, which is one of the C features that Frama-C's value analysis, as a static analyzer, does not claim to handle precisely in all cases.

Here, however, the program is completely deterministic: it does not have inputs of any kind. With an unreleased version of the value analysis, it is possible to interpret such programs for the detection of undefined behaviors. In this interpreter mode, the value analysis handles many more C constructs, including malloc() and realloc() calls.

That version of the value analysis warns about the use of p in the line if (p == q), which is indeed the first use of p after its contents have become indeterminate (in the sense of the C99 standard). The issue here relates to that post on the impossibility of freeing a pointer twice.

My entry

My entry was on the subject of illegal pointer arithmetic, which the value analysis does not warn about because there would be too many programs that it would (correctly) flag. It also involved comparisons of such illegal pointers, which the value analysis does warn about. These comparisons are also commonly found in the wild, but these are very very often actually dangerous, and have security implications.

There is nothing much I can say on the subject that was not already in a previous post on pointer comparisons.

Conclusion

You might think I am going to end up with a statement like “Frama-C's value analysis detects* all undefined behaviors because it finds the undefined behaviors in all four examples illustrated in the contest”.

Perish the thought! That would be fallacious. My claim is that Frama-C's value analysis detects* all undefined behaviors that you may care about, and even some you may think you do not need to care about but that might have surprising consequences, as exemplified on some notable entries in John's contest. The reason it finds* all undefined behaviors is not that it does well for a couple of examples but that it was designed to do so.

Incidentally, Frama-C's value analysis was tested on a lot more examples when working on this article on testcase reduction, and results were compared to KCC's. The Csmith random C program generator produces only defined programs, but C-Reduce, one of the reducers described in the article, produces tons of undefined behaviors. Although that was not a goal of the article, the work also served as an application of differential testing to the detectors of undefined behaviors that are KCC and Frama-C's value analysis.

There are problematic features of the C language that were left aside when initially designing the value analysis for its primary target of critical embedded code. This limits some uses, but it remains possible to extract valuable information off a C program, using the value analysis and the right methodology. One workaround is to use the value analysis as a C interpreter, in which case standard library functions are easier to model. Many of them already are available in the unreleased value analysis.

(*) when used properly. See Winner #1 and Winner #2 for examples of what “using properly” may entail.

My thanks to Sven Mattsen for proofreading this post.

I have been able to divert a few hours yesterday and today for programming. It was well worth it, as I have discovered a theorem. It is new to me, and I wonder whether it was ever published. The theorem is, a C program cannot double free() a block even if it tries! I thought this deserved an announcement. Some students I know would be glad to hear that.

Indeed, a C program can at most pass an indeterminate value to a function. This is highly objectionable in itself. Students, don't write programs that do this either!

Allow me to demonstrate on the following program:

...
main(){
  int *p = malloc(12);
  int *q = p;

  free(p);
  free(q);
}

You might think that the program above double-frees the block allocated and referenced by p. But it doesn't: after the first call to free(), the contents of q are indeterminate, so that it's impossible to free the block a second time. I will insert a few calls to primitive function Frama_C_dump_each() to show what happens:

main(){
  int *p = malloc(12);
  int *q = p;

  Frama_C_dump_each();
  free(p);
  Frama_C_dump_each();
  free(q);
  Frama_C_dump_each();
}

$ frama-c -val this_is_not_a_double_free.c

The analysis follows the control flow of the program. First a block is allocated, and a snapshot of the memory state is printed a first time in the log:

...
[value] computing for function malloc <- main.
        Called from test.c:12.
[value] Recording results for malloc
[value] Done for function malloc
[value] DUMPING STATE of file test.c line 15
        p ∈ {{ &Frama_C_alloc }}
        q ∈ {{ &Frama_C_alloc }}
        =END OF DUMP==
...

Then free() is called and a second snapshot is logged:

...
[value] computing for function free <- main.
        Called from test.c:16.
[value] Recording results for free
[value] Done for function free
[value] DUMPING STATE of file test.c line 17
        p ∈ ESCAPINGADDR
        q ∈ ESCAPINGADDR
        =END OF DUMP==
...

And then, at the critical moment, the program passes an indeterminate value to a function:

...
test.c:18:[kernel] warning: accessing left-value q that contains escaping addresses;
                      assert(Ook)
test.c:18:[kernel] warning: completely undefined value in {{ q -> {0} }} (size:<32>).
test.c:18:[value] Non-termination in evaluation of function call expression argument
        (void *)q
...
[value] Values at end of function main:
          NON TERMINATING FUNCTION

Note how the last call to Frama_C_dump_each() is not even printed. Execution stops at the time of calling the second free(). Generalizing this reasoning to all attempts to call free() twice, we demonstrate it is impossible to double-free a memory block.

QED

Question: does it show too much that this detection uses the same mechanism as the detection for addresses of local variables escaping their scope? I can change the message if it's too confusing. What would a better formulation be?

As I hope I have made clear, the message will be displayed as soon as the program does anything at all with a pointer that has been freed — and indeed, I have not implemented any special check for double frees. The program below is forbidden just the same, and does not get past the q++;

main(){
  int *p = malloc(12);
  int *q = p;

  Frama_C_dump_each();
  free(p);
  q++;
  ...
}

On the other hand, I still have to detect that the program does not free a global variable, something the analysis currently allows to great comical effect. (Hey, I have only had a couple of hours!)

Please leave your thoughts in the comments.

Introduction

In the first post in the obfuscated animated donut series, my colleague Anne pointed out that:

The alarm about : assert \valid(".,-~:;=!*#$@"+tmp_7); seems strange because the analysis tells us that tmp_7 ∈ [0..40] at this point... How can this be valid ?

It is only safe to use an offset between 0 and 11 when accessing inside the quantization string ".,-~:;=!*#$@". The value analysis determined that the offset was certainly between 0 and 40, so the analyzer can't be sure that the memory access is safe, and emits an alarm. Of course, [0..40] is only a conservative approximation for the value of the offset at that point. This approximation is computed with model functions sin() and cos() only known to return values in [-1.0 .. 1.0]. These model functions completely lose the relation between the sine and the cosine of a same angle. The computation we are interested in is more clearly seen in the function below, which was extracted using Frama-C's slicing plug-in.

/*@ ensures 0 <= \result <= 11 ; */
int compute(float A, float B, float i, float j)
{
  float c=sin(i), l=cos(i);
  float d=cos(j), f=sin(j);
  float g=cos(A), e=sin(A);
  float m=cos(B), n=sin(B);

  int N=8*((f*e-c*d*g)*m-c*d*e-f*g-l*d*n);
  return N>0?N:0;
}

Both the sine and the cosine of each argument are computed, and eventually used in the same formula. Surely the analyzer would do better if it could just remember that sin(x)²+cos(x)²=1. Alas, it doesn't. It just knows that each is in [-1.0 .. 1.0].

Making progress

But the Nitrogen release will have more precise sin() and cos() modelizations! Focusing only on the sub-computation we are interested in, we can split the variation domains for arguments A, B, i, j of function compute(), until it is clear whether or not N may be larger than 11.

I suggest to use the program below (download), where the code around compute() creates the right analysis context. Just like in unit testing, there are stubs and witnesses in addition to the function under test.

...

/*@ ensures 0 <= \result <= 11 ; */
int compute(float A, float B, float i, float j)
{
  float c=sin(i), l=cos(i);
  float d=cos(j), f=sin(j);
  float g=cos(A), e=sin(A);
  float m=cos(B), n=sin(B);

  int N=8*((f*e-c*d*g)*m-c*d*e-f*g-l*d*n);
Frama_C_dump_each();
  return N>0?N:0;
}

main(){
  int An, Bn, in, jn;
  float A, B, i, j;
  for (An=-26; An<26; An++)
    {
      A = Frama_C_float_interval(An / 8., (An + 1) / 8.);
      for (Bn=-26; Bn<26; Bn++)
	{
	  B = Frama_C_float_interval(Bn / 8., (Bn + 1) / 8.);
	  for (in=-26; in<26; in++)
	    {
	      i = Frama_C_float_interval(in / 8., (in + 1) / 8.);
	      for (jn=-26; jn<26; jn++)
		{
		  j = Frama_C_float_interval(jn / 8., (jn + 1) / 8.);
		  compute(A, B, i, j);
		}
	    }
	}
    }
}

Each for loop only covers the range [-3.25 .. 3.25]. In the real program, A and B vary in much wider intervals. If target functions sin() and cos() use a decent argument reduction algorithm (and incidentally, use the same one), the conclusions obtained for [-3.25 .. 3.25] should generalize to all finite floating-point numbers. This is really a property of the target platform that we cannot do anything about at this level, other than clearly state that we are assuming it.

In fact, the programs I analyzed were not exactly like above: this task is trivial to run in parallel, so I split the search space in four and ran four analyzers at the same time:

$ for i in 1 2 3 4 ; do ( time bin/toplevel.opt -val share/builtin.c \
an$i.c -obviously-terminates -no-val-show-progress \
-all-rounding-modes ) > log$i 2>&1 & disown %1 ; done

After a bit less than two hours (YMMV), I got files log1,..., log4 that contained snapshots of the memory state just after variable N had been computed:

[value] DUMPING STATE of file an1.c line 24
...
        c ∈ [-0.016591893509 .. 0.108195140958]
        l ∈ [-1. .. -0.994129657745]
        d ∈ [-1. .. -0.994129657745]
        f ∈ [-0.016591893509 .. 0.108195140958]
        g ∈ [-1. .. -0.994129657745]
        e ∈ [-0.016591893509 .. 0.108195140958]
        m ∈ [-1. .. -0.994129657745]
        n ∈ [-0.016591893509 .. 0.108195140958]
        N ∈ {-1; 0; 1}
        An ∈ {-26}
        Bn ∈ {-26}
        in ∈ {-26}
        jn ∈ {-26}
        A ∈ [-3.25 .. -3.125]
        B ∈ [-3.25 .. -3.125]
        i ∈ [-3.25 .. -3.125]
        j ∈ [-3.25 .. -3.125]
        =END OF DUMP==
an1.c:15:[value] Function compute: postcondition got status valid.
[value] DUMPING STATE of file an1.c line 24
...
        c ∈ [-0.016591893509 .. 0.108195140958]
        l ∈ [-1. .. -0.994129657745]
        d ∈ [-0.999862372875 .. -0.989992439747]
        f ∈ [-0.141120016575 .. -0.0165918916464]
        g ∈ [-1. .. -0.994129657745]
        e ∈ [-0.016591893509 .. 0.108195140958]
        m ∈ [-1. .. -0.994129657745]
        n ∈ [-0.016591893509 .. 0.108195140958]
        N ∈ {-2; -1; 0; 1}
        An ∈ {-26}
        Bn ∈ {-26}
        in ∈ {-26}
        jn ∈ {-25}
        A ∈ [-3.25 .. -3.125]
        B ∈ [-3.25 .. -3.125]
        i ∈ [-3.25 .. -3.125]
        j ∈ [-3.125 .. -3.]
        =END OF DUMP==
[value] DUMPING STATE of file an1.c line 24
...

Some input sub-intervals for A, B, i, j make the analyzer uncertain about the post-condition N ≤ 11. Here are the first such inputs, obtained by searching for "unknown" in the log:

...
[value] DUMPING STATE of file an1.c line 24
 ...
        c ∈ [-0.850319802761 .. -0.778073191643]
        l ∈ [-0.628173649311 .. -0.526266276836]
        d ∈ [0.731688857079 .. 0.810963153839]
        f ∈ [0.585097253323 .. 0.681638777256]
        g ∈ [-1. .. -0.994129657745]
        e ∈ [-0.016591893509 .. 0.108195140958]
        m ∈ [-1. .. -0.994129657745]
        n ∈ [-0.016591893509 .. 0.108195140958]
        N ∈ {8; 9; 10; 11; 12}
        An ∈ {-26}
        Bn ∈ {-26}
        in ∈ {-18}
        jn ∈ {5}
        A ∈ [-3.25 .. -3.125]
        B ∈ [-3.25 .. -3.125]
        i ∈ [-2.25 .. -2.125]
        j ∈ [0.625 .. 0.75]
        =END OF DUMP==
an1.c:15:[value] Function compute: postcondition got status unknown.
...

Here is the complete list of value sets obtained for variable N on the entire search space:

$ cat log? | grep "N " | sort | uniq
        N ∈ {0; 1}
        N ∈ {0; 1; 2}
        N ∈ {0; 1; 2; 3}
        N ∈ {0; 1; 2; 3; 4}
        N ∈ {0; 1; 2; 3; 4; 5}
        N ∈ {0; 1; 2; 3; 4; 5; 6}
        N ∈ {-1; 0}
        N ∈ {-1; 0; 1}
        N ∈ {10; 11}
        N ∈ {10; 11; 12}
        N ∈ {-1; 0; 1; 2}
        N ∈ {-1; 0; 1; 2; 3}
        N ∈ {-1; 0; 1; 2; 3; 4}
        N ∈ {-10; -9}
        N ∈ {-10; -9; -8}
        N ∈ {-10; -9; -8; -7}
        N ∈ {-10; -9; -8; -7; -6}
        N ∈ {-10; -9; -8; -7; -6; -5}
        N ∈ {-10; -9; -8; -7; -6; -5; -4}
        N ∈ {-11; -10}
        N ∈ {-11; -10; -9}
        N ∈ {-11; -10; -9; -8}
        N ∈ {-11; -10; -9; -8; -7}
        N ∈ {-11; -10; -9; -8; -7; -6}
        N ∈ {-11; -10; -9; -8; -7; -6; -5}
        N ∈ {1; 2}
        N ∈ {-12; -11; -10}
        N ∈ {-12; -11; -10; -9}
        N ∈ {-12; -11; -10; -9; -8}
        N ∈ {-12; -11; -10; -9; -8; -7}
        N ∈ {-12; -11; -10; -9; -8; -7; -6}
        N ∈ {1; 2; 3}
        N ∈ {1; 2; 3; 4}
        N ∈ {1; 2; 3; 4; 5}
        N ∈ {1; 2; 3; 4; 5; 6}
        N ∈ {1; 2; 3; 4; 5; 6; 7}
        N ∈ {-13; -12; -11; -10; -9}
        N ∈ {-13; -12; -11; -10; -9; -8}
        N ∈ {-13; -12; -11; -10; -9; -8; -7}
        N ∈ {-2; -1}
        N ∈ {-2; -1; 0}
        N ∈ {-2; -1; 0; 1}
        N ∈ {-2; -1; 0; 1; 2}
        N ∈ {-2; -1; 0; 1; 2; 3}
        N ∈ {2; 3}
        N ∈ {2; 3; 4}
        N ∈ {2; 3; 4; 5}
        N ∈ {2; 3; 4; 5; 6}
        N ∈ {2; 3; 4; 5; 6; 7}
        N ∈ {2; 3; 4; 5; 6; 7; 8}
        N ∈ {-3; -2}
        N ∈ {-3; -2; -1}
        N ∈ {-3; -2; -1; 0}
        N ∈ {-3; -2; -1; 0; 1}
        N ∈ {-3; -2; -1; 0; 1; 2}
        N ∈ {3; 4}
        N ∈ {3; 4; 5}
        N ∈ {3; 4; 5; 6}
        N ∈ {3; 4; 5; 6; 7}
        N ∈ {3; 4; 5; 6; 7; 8}
        N ∈ {3; 4; 5; 6; 7; 8; 9}
        N ∈ {-4; -3}
        N ∈ {-4; -3; -2}
        N ∈ {-4; -3; -2; -1}
        N ∈ {-4; -3; -2; -1; 0}
        N ∈ {-4; -3; -2; -1; 0; 1}
        N ∈ {4; 5}
        N ∈ {4; 5; 6}
        N ∈ {4; 5; 6; 7}
        N ∈ {4; 5; 6; 7; 8}
        N ∈ {4; 5; 6; 7; 8; 9}
        N ∈ {4; 5; 6; 7; 8; 9; 10}
        N ∈ {-5; -4}
        N ∈ {-5; -4; -3}
        N ∈ {-5; -4; -3; -2}
        N ∈ {-5; -4; -3; -2; -1}
        N ∈ {-5; -4; -3; -2; -1; 0}
        N ∈ {5; 6}
        N ∈ {5; 6; 7}
        N ∈ {5; 6; 7; 8}
        N ∈ {5; 6; 7; 8; 9}
        N ∈ {5; 6; 7; 8; 9; 10}
        N ∈ {5; 6; 7; 8; 9; 10; 11}
        N ∈ {-6; -5}
        N ∈ {-6; -5; -4}
        N ∈ {-6; -5; -4; -3}
        N ∈ {-6; -5; -4; -3; -2}
        N ∈ {-6; -5; -4; -3; -2; -1}
        N ∈ {-6; -5; -4; -3; -2; -1; 0}
        N ∈ {6; 7}
        N ∈ {6; 7; 8}
        N ∈ {6; 7; 8; 9}
        N ∈ {6; 7; 8; 9; 10}
        N ∈ {6; 7; 8; 9; 10; 11}
        N ∈ {6; 7; 8; 9; 10; 11; 12}
        N ∈ {-7}
        N ∈ {7}
        N ∈ {-7; -6}
        N ∈ {-7; -6; -5}
        N ∈ {-7; -6; -5; -4}
        N ∈ {-7; -6; -5; -4; -3}
        N ∈ {-7; -6; -5; -4; -3; -2}
        N ∈ {-7; -6; -5; -4; -3; -2; -1}
        N ∈ {7; 8}
        N ∈ {7; 8; 9}
        N ∈ {7; 8; 9; 10}
        N ∈ {7; 8; 9; 10; 11}
        N ∈ {7; 8; 9; 10; 11; 12}
        N ∈ {7; 8; 9; 10; 11; 12; 13}
        N ∈ {-8; -7}
        N ∈ {-8; -7; -6}
        N ∈ {-8; -7; -6; -5}
        N ∈ {-8; -7; -6; -5; -4}
        N ∈ {-8; -7; -6; -5; -4; -3}
        N ∈ {-8; -7; -6; -5; -4; -3; -2}
        N ∈ {8; 9}
        N ∈ {8; 9; 10}
        N ∈ {8; 9; 10; 11}
        N ∈ {8; 9; 10; 11; 12}
        N ∈ {8; 9; 10; 11; 12; 13}
        N ∈ {9; 10}
        N ∈ {9; 10; 11}
        N ∈ {9; 10; 11; 12}
        N ∈ {9; 10; 11; 12; 13}
        N ∈ {-9; -8}
        N ∈ {-9; -8; -7}
        N ∈ {-9; -8; -7; -6}
        N ∈ {-9; -8; -7; -6; -5}
        N ∈ {-9; -8; -7; -6; -5; -4}
        N ∈ {-9; -8; -7; -6; -5; -4; -3}

The above is pretty cute, and it is just short enough for manual inspection. The inspection reveals that all the lines with dangerous values for N contain the number 12, so that we can grep for that. We arrive at the following list of dangerous value sets for N:

        N ∈ {10; 11; 12}
        N ∈ {6; 7; 8; 9; 10; 11; 12}
        N ∈ {7; 8; 9; 10; 11; 12}
        N ∈ {7; 8; 9; 10; 11; 12; 13}
        N ∈ {8; 9; 10; 11; 12}
        N ∈ {8; 9; 10; 11; 12; 13}
        N ∈ {9; 10; 11; 12}
        N ∈ {9; 10; 11; 12; 13}

Note: it is normal that each dangerous line contains 12: if the safety property that we are trying to establish is true, then there must be a safe value (say, 10) in each set. And the sets contain several values as the result of approximations in interval arithmetics, so that sets of values for N are intervals of integers. Any interval that contains both a safe value and an unsafe value must contain the largest safe value 11 and the smallest unsafe value 12.

Conclusion

We have improved our estimate of the maximum value of N from 40 to 13. That's pretty good, although it is not quite enough to conclude that the donut program is safe.

In conclusion to this post, now that we have ascertained that a dangerous set for N was a set that contains 12, we can compute the ratio of the search space for which we haven't yet established that the program was safe. We are exploring a dimension-4 cube by dividing it in sub-cubes. There are 52⁴ = 7 311 616 sub-cubes. Of these, this many are still not known to be safe:

$ cat log? | grep "N " | grep " 12" | wc -l
360684

We are 95% there. The remaining 5% shouldn't take more than another 95% of the time.

This post is about widening. This technique was shown in the second part of a previous post about memcpy(), where it was laboriously used to analyze imprecisely function memcpy() as it is usually written. The value analysis in Frama-C has the ability to summarize loops in less time than they would take to execute. Indeed, it can analyze in finite time loops for which it is not clear whether they terminate at all. This is the result of widening, the voluntary loss of precision in loops. As a side-effect, widening can produce some counter-intuitive results. This post is about these counter-intuitive results.

The expected: better information in, better information out

Users generally expect that the more you know about the memory state before the interpretation of a piece of program, the more you should know about the memory state after this piece of program. Users are right to expect this, and indeed, it is generally true:

int u;

main(){
  u = Frama_C_interval(20, 80);
  if (u & 1)
    u = 3 * u + 1; // odd
  else
    u = u / 2; // even
}

The command frama-c -val c.c share/builtin.c shows the final value of u to be in [10..241]. If the first line of main() is changed to make the input value range over [40..60] instead of [20..80], then the output interval for u becomes [20..181]. This seems normal : you know more about variable u going in (all the values in [40..60] are included in [20..80]), so you know more about u coming out ([20..181] is included in [10..241]). Most of the analyzer's weaknesses, such as, in this case, the inability to recognize that the largest value that can get multiplied by 3 is respectively 79 or 59, do not compromise this general principle. But widening does.

Better information in, worse information out

Widening happens when there is a loop in the program, such as in the example below:

int u, i;

main(){  
  u = Frama_C_interval(40, 60);
  for (i=0; i<10; i++)
    u = (25 + 3 * u) / 4;
}

$ frama-c -val loop.c share/builtin.c
...
[value] Values for function main:
          u ∈ [15..60]

Analyzing this program with the default options produce the interval [15..60] for u. Now, if we change the initial set for u to be [25..60], an interval that contains all the values from [40..60] and more:

$ frama-c -val loop.c share/builtin.c
...
[value] Values for function main:
          u ∈ [25..60]

By using a larger, less precise set as the input of this analysis, we obtain a more precise result. This is counter-intuitive, although it is rarely noticed in practice.

Better modelization, worse analysis

There are several ways to make this unpleasant situation less likely to be noticed. Only a few of these are implemented in Frama-C's value analysis. Still, it was difficult enough to find a short example to illustrate my next point, and the program is only as short as I was able to make it. Consider the following program:

int i, j;
double ftmp1;

double Frama_C_cos(double);
double Frama_C_cos_precise(double);

#define TH 0.196349540849361875
#define NBC1 14
#define S 0.35355339

double cos(double x)
{
  return Frama_C_cos(x);
}

main(){  
  for (i = 0; i < 8; i++)
    for (j = 0; j < 8; j++)
      {
        ftmp1 = ((j == 0) ? S : 0.5) *
          cos ((2.0 * i + 1.0) * j * TH);
	ftmp1 *= (1 << NBC1);
      }
}

If you analyze the above program with a smartish modelization of function cos(), one that recognizes when its input is an interval [l..u] of floating-point values on which cos() is monotonous, and computes the optimal output in this case, you get the following results:

$ frama-c -val idct.c share/builtin.c
...
          ftmp1 ∈ [-3.40282346639e+38 .. 8192.]

Replacing the cos() modelization by a more naive one that returns [-1. .. 1.] as soon as its input set isn't a singleton, you get the improved result:

$ frama-c -val idct.c share/builtin.c
...
          ftmp1 ∈ [-8192. .. 8192.]

This is slightly different from the [25..60] example, in that here we are not changing the input set but the modelization of one of the functions involved in the program. Still, the causes are the same, and the effects just as puzzling when they are noticeable. And this is why, until version Carbon, only the naive modelization of cos() was provided.

In the next release, there will be, in addition to the naive versions, Frama_C_cos_precise() and Frama_C_sin_precise() built-ins for when you want them. Typically, this will be for the more precise analyses that do not involve widening.

This post follows that post where a typical double-precision implementation for a cosine function is shown. That implementation is not a correctly rounded implementation, but its double result is precise to a level close to what the double-precision floating-point representation allows. It uses low-level tricks that make it faster on nineties computers.

In this post, I offer a slightly simplified version of the same function (download). This function does not compute exactly the same results, but according to a few billion tests I ran, it does not seem more wrong than the original. I expect the simplified version is actually faster than the original on tens computers, but this is not why I modified it. I needed to make changes because the value analysis does not yet precisely handle all the low-level conversions in the original (accessing the 32 most significant bits of a double). Although the value analysis does not say anything false, it is only possible to analyze the function precisely if this idiom is handled precisely.

If it was important to analyze as-is code that uses this low-level idiom, implementation of this treatment could be prioritized according to the results of the experiment in this post. Results on the original function can be expected to be the same after the new feature is implemented as the results for the modified function described here.

I implemented a quick-and-dirty value analysis builtin for the purpose of this verification. It is not provided so much for being re-used as for possible comments: it is large enough that a bug in it could invalidate the verification. This builtin checks that its second argument is close enough to the cosine of its first argument (it uses its third argument as allowable distance). The builtin displays a log line for each call, so that you can check that it was properly called.

Demonstrating how easy it is to write an unsafe builtin, this implementation assumes without verification that it is called in a range where function cosine is decreasing.

I then wrote the following lines of C:

          x = Frama_C_double_interval(..., ...);
          m = mcos(x);
          Frama_C_compare_cos(x, m, 0x1.0p-21);

Next, I arranged to have the value analysis repeatedly analyze these lines with bounds passed to Frama_C_double_interval that, taken together, cover entirely the range from 0 to π/4. I will not say how I did, because it is not very interesting. Besides, the method I used is convenient enough to discourage someone from writing an Emacs Lisp macro to generate sub-division assertions. That would be an unwanted side-effect: it in fact has complementary advantages with the would-be assertion generation macro that I recommend someone look into.

The verification used a few hours of CPU time (I did not try very hard to make it fast). One thing I did was adapt the input sub-interval width to the part of the curve that was being verified, from two millionths for inputs near zero to a quarter of a millionth for inputs near π/4.

During this experiment, I actually used the cosine implementation from CRlibm as reference function. This library, incidentally developed at my alma mater, provides a number of fast and correctly rounded implementations for mathematical functions. The value analysis was configured with the -all-rounding-modes option, so the results are guaranteed for arbitrary rounding modes, and are guaranteed even if the compiler uses supplemental precision for intermediate results (e.g. if the compiler uses the infamous 80-bit historial floating-point unit of Intel processors).

Here are the logs I obtained (warning: I did not make any file of that level of boringness publicly available since my PhD, plus these are really big, even after compression. Handle with care): log1.gz log2.gz log3.gz log4.gz log5.gz log6.gz log7.gz log8.gz log9.gz log10.gz log11.gz log12.gz log13.gz log14.gz log15.gz. To save you the bother of downloading even the first file, here is what the first three lines look like.

[value] CC -0.0000000000000000 0.0000019073486328 0.9999999999981809 1.0000000000000002 0.9999999999981810 1.0000000000000000 OK
[value] CC 0.0000019073486328 0.0000038146972656 0.9999999999927240 0.9999999999981811 0.9999999999927240 0.9999999999981811 OK
[value] CC 0.0000038146972656 0.0000057220458984 0.9999999999836291 0.9999999999927242 0.9999999999836291 0.9999999999927242 OK

From the fact that the builtin did not complain for any sub-interval, we can conclude that on the interval 0..π/4, the proposed implementation for cosine is never farther than 2^-21 (half a millionth) from the real function cosine.

According to the aforementioned billion tests I made, the implementation seems to be much more precise than above (to about 2^-52 times its argument). The key difference is that I am quite sure it is precise to half a millionth, and that I have observed that it seems to be precise to about 2^-52 times its argument.

Guillaume Melquiond provided irreplaceable floating-point expertise in the preparation of this post.