Using NativeCall to call the C fn `erf` gets more precise output than `erf` in C

174 Views Asked by At

I have written a Raku script to call erf function in C standard library:

use NativeCall;
sub erf(num64) returns num64 is native { * };

say [0.5,1,2,3,4,-0.9].map: {erf($_.Num)};

The output of this script

(0.5204998778130465 0.8427007929497149 0.9953222650189527 0.9999779095030014 0.9999999845827421 -0.7969082124228322)

matches with the output from C for all values [0.5,1,2,3,4,-0.9] except 4.

For 4 the C outputs 1.000000 while Raku gives 0.9999999845827421.

To test the output for 4 in C, run this code:

#include <stdio.h> // Including header file for printf function
#include <math.h>  // Including header file for erf function

int main (){

  double param, result;
  param = 4.0;
  result = erf(param);
  printf("erf (%f) = %f\n", param, result);
  return 0;
}

Any idea what's going on? I need to output 1.0 from Raku too.

3

There are 3 best solutions below

0
Elizabeth Mattijsen On BEST ANSWER

You're comparing apples with oranges.

use NativeCall;
sub erf(num64) returns num64 is native { * };

say .fmt("%f") for [0.5,1,2,3,4,-0.9].map: {erf($_.Num)}

0.520500
0.842701
0.995322
0.999978
1.000000
-0.796908

You're using printf in C, if you use .fmt (the easier way to say sprintf in Raku), then you'd also get 1.0.

0
Eric Postpischil On

For the C code, change %f to %.99g to show more digits. This reveals erf(4) returns 0.9999999845827420852373279558378271758556365966796875.

%f requests six digits after the decimal point. The value is rounded to fit that format. %.numberf requests number digits after the decimal point and always used the “fixed” format. %.numberg requests number significant digits and uses a a “general” format that switches to exponential notation when appropriate.

For the Raku code, if you want output of “1.0” or “1.000000”, you will need to apply some formatting request to the output. I do not practice Raku, but a brief search shows Raku has printf-like features you can use, so requesting the %f format with it should duplicate the C output.

2
Ted Lyngmo On

Yet another version building on @Eric's explanation about precision:

my $fmt = '%.6f';
printf [$fmt ~ ' '] x 5 ~ $fmt ~ "\n", [0.5,1,2,3,4,-0.9].map: {erf($_.Num)};

[$fmt ~ ' '] x 5 ~ $fmt ~ "\n" builds the format string "%.6f %.6f %.6f %.6f %.6f %.6f\n"

Output

0.520500 0.842701 0.995322 0.999978 1.000000 -0.796908

Or even more as proposed by raiph, using the infix operator xx:

printf "{'%.6f' xx 6} \n", [0.5,1,2,3,4,-0.9].map: {erf($_.Num)};