summaryrefslogtreecommitdiff
diff options
context:
space:
mode:
-rw-r--r--fourier.c90
1 files changed, 36 insertions, 54 deletions
diff --git a/fourier.c b/fourier.c
index af7c5e8..a054c8f 100644
--- a/fourier.c
+++ b/fourier.c
@@ -35,9 +35,6 @@ static inline double TrapezoidIntegrate(double x0,
double x1,
int n,
function_t select);
-static inline double thefunction(double x,
- int n,
- function_t select);
/**************
** DoFourier **
@@ -72,7 +69,7 @@ DoFourier(void)
** less than or equal to the permitted minimum, re-allocate
** larger arrays and try again.
*/
- } while (DoFPUTransIteration(abase,bbase, array_size) <= MINIMUM_TICKS);
+ } while (DoFPUTransIteration(abase, bbase, array_size) <= MINIMUM_TICKS);
} else {
/*
** Don't need self-adjustment. Just allocate the
@@ -99,8 +96,8 @@ DoFourier(void)
*************************
** Perform an iteration of the FPU Transcendental/trigonometric
** benchmark. Here, an iteration consists of calculating the
-** first n fourier coefficients of the function (x+1)^x on
-** the interval 0,2. n is given by arraysize.
+** first NUM_STEPS fourier coefficients of the function (x+1)^x on
+** the interval 0,2.
*/
static clock_t
DoFPUTransIteration(double *abase, double *bbase, unsigned long arraysize)
@@ -151,56 +148,41 @@ static inline double
TrapezoidIntegrate(double x0, double x1, int n, function_t select)
{
double dx = (x1 - x0) / (double)NUM_STEPS; /* Stepsize */
- double rvalue;
- int num_steps = NUM_STEPS - 1; /* Already done 1 step */
+ double rvalue = pow(x0 + 1.0, x0);
+ int num_steps = NUM_STEPS - 2; /* Already done 1 step */
- rvalue = thefunction(x0, n, select);
- /*
- * Compute the other terms of the integral.
- */
- while(--num_steps) {
- x0 += dx;
- rvalue += thefunction(x0, n, select);
- }
+ switch (select) {
+ case NONE:
+ while(num_steps) {
+ x0 += dx;
+ rvalue += pow(x0 + 1.0, x0);
+ num_steps--;
+ }
- /*
- * Finish computation
- */
- rvalue += thefunction(x1, n, select);
+ rvalue += pow(x1 + 1.0, x1);
+ break;
+ case COS:
+ rvalue *= cos(M_PI * n * x0);
+
+ while(num_steps) {
+ x0 += dx;
+ rvalue += pow(x0 + 1.0, x0) * cos(M_PI * n * x0);
+ num_steps--;
+ }
+ rvalue += pow(x1 + 1.0, x1) * cos(M_PI * n * x1);
+ break;
+ case SIN:
+ rvalue *= sin(M_PI * n * x0);
+
+ while(num_steps) {
+ x0 += dx;
+ rvalue += pow(x0 + 1.0, x0) * sin(M_PI * n * x0);
+ num_steps--;
+ }
+
+ rvalue += pow(x1 + 1.0, x1) * sin(M_PI * n * x1);
+ break;
+ }
return rvalue / 2.0 * dx;
}
-
-/****************
-** thefunction **
-*****************
-** This routine selects the function to be used
-** in the Trapezoid integration.
-** x is the independent variable
-** n is the series number
-** select chooses which of the sine/cosine functions
-** are used. note the special case for select=0.
-*/
-static inline double thefunction(double x, /* Independent variable */
- int n, /* Omega * term */
- function_t select) /* Choose term */
-{
-
-/*
-** Use select to pick which function we call.
-*/
-switch(select)
-{
- case 0: return(pow(x+(double)1.0,x));
-
- case 1: return(pow(x+(double)1.0,x) * cos(M_PI * n * x));
-
- case 2: return(pow(x+(double)1.0,x) * sin(M_PI * n * x));
-}
-
-/*
-** We should never reach this point, but the following
-** keeps compilers from issuing a warning message.
-*/
-return(0.0);
-}