Originally Posted by

**aghast**
1. Your factorial is taking int parameters. You should print out the factorials of all the numbers up to 201! to be sure it works. (Hint: try it on your favorite calculator first.)

A small program to test this:

Code:

#include <stdio.h>
#include <math.h>
float fact ( unsigned int x )
{
float n = 1.0;
while ( x > 1 )
n *= x--;
return n;
}
int main ( void )
{
float n, m;
int i;
n = 1.0;
for ( i = 2; i <= 201; i++ )
{
m = fact ( i );
printf ( "%d -> %g ", i, m );
// stops if is a nan or inf.
if ( ! isnormal ( m ) )
{
puts ( "[fail]." );
break;
}
n = m;
putchar ( '\n' );
}
}

Running it you can see it fails for factorial of 35. Changing to double and will fail for fatorial of 171. It works for factorial of 201 with long double (will fail for factorial of 1755).

2. The expression -1 ^n is really just -1 * -1 * -1 ..., out to n times. Thus, -1^0 is 1, -1^1 is -1, -1^2 is 1, etc. These are just +1/-1 alternating each time.

As such it can be easily coded as:

Code:

int signal_( unsigned int x ) { return x & 1 ? -1 : 1; }

5. Since the factorials you are being asked to compute are increasing, you can get much better performance by remembering the last argument and result in your factorial function.

Or, since there are a fixed interations for the loop, use a pre-compiled table:

Code:

/* gentbl.c */
#include <stdio.h>
static long double fact ( unsigned int x )
{
long double n = 1.0;
while ( x > 1 )
n *= x--;
return n;
}
int main ( void )
{
unsigned int i;
fputs ( "long double tbl_[] = {\n", stdout );
for ( i = 0; i < 100; i++ )
printf ( " %.18LgL,\n", fact ( 2*i + 1 ) );
printf ( " %.18LgL\n};\n", fact ( 2*i + 1 ) );
}

Compile, and running to create the table:

Code:

$ cc -o gentable gentable.c
$ ./gentable > tbl.c

The sin_() function:

Code:

#include <stdio.h>
#include <math.h>
/* from tbl.c */
extern long double tbl_[];
double sin_( double x )
{
unsigned int i;
long double n;
double t;
// x cannot be a non-number
if ( isnan( x ) || isinf( x ) )
return NAN;
// Restrict x to (-2*PI,2*PI).
// not sure if this will have an impact
// on "precision".
t = fabs( x );
if ( t >= 2.0 * M_PI )
x /= t / (2.0 * M_PI );
for ( i = 0; i <= 100; i++ )
{
n += (i & 1 ? -1 : 1 ) * powl( x, 2*i + 1 ) / tbl_[i];
// powl() can result in an overflow!
if ( isinf( n ) || isnan( n ) )
break;
}
// n between -1 and 1 will fit in a double.
return n;
}