Octave & C

A bit of a complicated build process, but kind of a novel approach and the results were encouraging.

The approach is to generate approximating quadratic polynomial for each degree. So degree = [0, 1), degree = [1, 2), etc. each have a different polynomial.

Octave - building part

Octave is publicly available - Google download octave.

This determines the best fit quadratic curve for every degree.

Save as build-fast-trig.m:

format long;
for d = 0:359
    x = (d-1):0.1:(d+1);
    y = sin(x / 360 * 2 * pi);
    polyfit(x, y, 2)
endfor

C - building part

This converts doubles in text format to native binary format on your system.

Save as build-fast-trig.c:

#include <stdio.h>

int main()
{
    double d[3];

    while (scanf("%lf %lf %lf", d, d + 1, d + 2) == 3)
        fwrite(d, sizeof(double), 3, stdout);

    return 0;
}

Compile:

gcc -o build-fast-trig build-fast-trig.c

Generating the coefficients file

Run:

octave build-fast-trig.m | grep '^ ' | ./build-fast-trig > qcoeffs.dat

Now we have qcoeffs.dat as the data file to use for the actual program.

C - fast-trig part

Save as fast-trig.c:

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

#define INPUT    "qcoeffs.dat"

#define DEGREES    360

typedef struct {double a, b, c;} QCOEFFS;

double normalize(double d)
{
    if (d < 0.0)
        d += ceil(d / -(double)DEGREES) * (double)DEGREES;

    if (d >= DEGREES)
        d -= floor(d / (double)DEGREES) * (double)DEGREES;

    return d;
}

int main()
{
    FILE *f;
    time_t tm;
    double d;
    QCOEFFS qc[DEGREES];

    if (!(f = fopen(INPUT, "rb")) || fread(qc, sizeof(QCOEFFS), DEGREES, f) < DEGREES)
    {
        fprintf(stderr, "Problem with %s - aborting.", INPUT);
        return EXIT_FAILURE;
    }
    fclose(f);

    tm = -clock();

    while (scanf("%lf", &d) > 0)
    {
        int i;
        double e, f;

        /* sin */
        d = normalize(d);
        i = (int)d;
        e = (qc[i].a * d + qc[i].b) * d + qc[i].c;

        /* cos */
        d = normalize((double)DEGREES / 4.0 - d);
        i = (int)d;
        f = (qc[i].a * d + qc[i].b) * d + qc[i].c;

        /* tan = sin / cos */

        /* output - format closest to specs, AFAICT */
        if (d != 0.0 && d != 180.0)
            printf("%.6e %.6e %.6e\n", e, f, e / f);
        else
            printf("%.6e %.6e n\n", e, f);
    }

    tm += clock();

    fprintf(stderr, "time: %.3fs\n", (double)tm/(double)CLOCKS_PER_SEC);    

    return EXIT_SUCCESS;
}

Compile:

gcc -o fast-trig fast-trig.c -lm

Run:

./fast-trig < trig.in > trig.out

It will read from trig.in, save to trig.out and print to console the elapsed time with millisecond precision.

Depending on the testing methods used it may fail on certain input, e.g.:

$ ./fast-trig 
0
-6.194924e-19 1.000000e+00 -6.194924e-19

The correct output should be 0.000000e+00 1.000000e+00 0.000000e+00. If the results are validated using strings, the input will fail, if they are validated using absolute error, e.g. fabs(actual - result) < 1e-06, the input will pass.

The maximum absolute error for sin and cos was < 1e-08. For tan, because the result isn't limited to ± 1, the absolute error could be larger, but the maximum relative error, e.g. fabs((actual - result) / actual) for |actual| > 1 was < 2e-08.

In testing, the average time per single input was (1.053 ± 0.007) µs, which on my machine was about 0.070 µs faster than native sin and cos.