...And then, of course, I discover the program does not calculate f**(2*S) but f**(2*S-1). Which makes my "correct" results not so very correct. Ouch.
It outputs: nearest = 0.36787944637187158792 down = 0.36787944526181193261 while we have
correct = 0.36787944576993540330 so both rounding modes are about equal.
Now you can say that this is because x already is rounded from the true value,
... and this is indeed ...
the case: taking S=2**25, we get nearest = 0.36787944391225085860 down = 0.36787944204965455918 with
correct = 0.36787944391235777181 so in that case round towards zero gets twice as many bits wrong as round to nearest, just as a naive analysis would show. Segher