Question

I was wondering if increasing the machine precision from float to double in C automatically raise terms in series expansion e.g. trigonometric functions ?

This and my other questions on stack overflow has to do with investigating occurring numerical instability on a wave equation in anisotropic media.

Was it helpful?

Solution

Generally, what affects the accuracy of routines to perform sine, logarithm, and so on is which version of the routine you call. A good math library has separate routines for single-precision sine (C’s sinf function), double-precision sine (C’s sin), and long-double sine (C’s sinl). In C, you usually call these versions explicitly, by writing calls to sinf, sin, or sinl. C also offers <tgmath.h>, which causes the source code sin(x) to expand to a specific version depending on the type of x. In C++, a function call will also be resolved depending on the type of the variable.

In a good math library, the sinf routine will use a faster algorithm with accuracy suited to the precision of float, while sin will use a slower algorithm suited to the precision of double. The quality of math libraries varies, as writing these routines is a complicated task.

Series expansion is not used. (In particular, Taylor series are not used due to poor error distribution and requiring too many terms to converge.) Instead, carefully prepared approximating polynomials are used. Some form of minimax polynomial is often used. A routine for a more precise type is likely to use a polynomial with more terms, but it is also likely to change in other ways, such as partitioning the domain into more intervals or using some form of extended precision. None of this is automatic; the routines are prepared manually by software enginers.

OTHER TIPS

According to this reference, while there are different versions of the sin(x) function for C++,

In C, only the double version of this function exists with this name

To confirm this, I wrote the following lines of code:

#include <stdio.h>
#include <math.h>

int main(void){

    double d = 0.12345e16;
    float f;
    f = d;

    printf("the difference is %f\n", sin(d) - sin(f));

    f = 0.12345e16;
    d = f;
    printf("the difference is now %f\n", sin(d) - sin(f));

}

This produces the following output (when compiled with the C compiler):

the difference is 1.947785
the difference is now 0.000000

EDIT - had a typo in my original code. Now updated, and giving the expected result.

The above shows that the two arguments are evaluated the same way, confirming that the same algorithm is used.

This is using the gcc compiler on Linux, with the -lm option, and with the file called fltchk.c (which I believe should cause gcc to use the C compiler).

Licensed under: CC-BY-SA with attribution
Not affiliated with StackOverflow
scroll top