#include <stdio.h>      /* standard I/O                */
#include <string.h>     /* for strcpy - 3 occurrences  */
#include <stdlib.h>     /* for exit   - 1 occurrence   */
#include <math.h>       /* for sin, exp etc.           */
/* #define enable_threading 1 */   /* uncomment for threading */
#ifdef enable_threading
#include <pthread.h>
#endif

struct timeval starttime,endtime;
double te0,te1,te2,te3,te4,te5,te6,te7;
double e1[4096]; /* 64MB array with 8 byte double systems */
long ae = 4096;
int n;
double ip;
double t = 0.49999975;
double t0 = 0.49999975;
double t1 = 0.50000025;
double t2 = 2.0;
double t3 = 1297263.1415998;
double lsfpu;
double x,y,z;
long i,ix,j,k,l;
double pass0,pass1,pass2,pass3,pass4,pass5,pass6,pass7;
long x100 = 222;
long xtra = 100;
long n1;
long n2;
long n3;
long n4;
long n5;
long n6;
long n7;
long n8;
long n1mult;
void pa(double e[4], double t, double t2);	
void p3(double *x, double *y, double *z, double t, double t1, double t2);
void po(double e1[4], long j, long k, long l);

double Section1() {
	        /* Section 1, Math: general math.h capabilities*/
	/* double_size=sizeof(double); */
        /* asize=1024*1024/double*mt; how many longs then in one array? */
	printf("\n\tPreparing test. Setting up a 64MB array and assigning arbitrary values....");
	for (ix=0; ix<(ae-1); ix++) {
		e1[ix] = 1.3333;
		e1[ix+1] = 8.9273;
		e1[ix+2] = 2361.24521;
		e1[ix+3] = 902521468.102937;
		e1[ix+4] = 5823.443;
                e1[ix+5] = 17.8;
                e1[ix+6] = 281621463.4121727;
                e1[ix+7] = 77321.81;
		ix=ix + 7;
	}
	printf("\n\tdone.");
	printf("\n\tTest running...\n");
	                /* Simple add/sub/mult/div functions 511 passes to fill array 8 at a time to increase randomness of numbers*/
       	gettimeofday(&starttime, NULL);
	for (i=0; i<x100; i++) {
#ifdef enable_threading
#pragma omp parallel for
#endif
        	for (ix=0; ix<(ae-1); ix++) 
                {
                e1[ix] = e1[ix] + e1[ix] - e1[ix] * e1[ix] / e1[ix] + e1[ix] - e1[ix] * e1[ix] / e1[ix] + t0 - t1 / t2 * t3;
                }
	}
	gettimeofday(&endtime, NULL);
        te0=((double)(endtime.tv_sec*1000000-starttime.tv_sec*1000000+endtime.tv_usec-starttime.tv_usec))/1000000;	
	te1 = 13 * ae * x100 / te0;
        printf("\n\tTest #1:\t Math.h add/sub/mult/div:\t%lf/s", te1);
		/* Trig functions */
	gettimeofday(&starttime, NULL);
	for (i=0; i<x100; i++) {
#ifdef enable_threading
#pragma omp parallel for
#endif
		for (ix=0; ix<(ae-1); ix++)
              {
                 e1[ix] = cos(e1[ix]) + sin(e1[ix]) - tan(e1[ix]) * acos(e1[ix]) / asin(e1[ix]) + atan(e1[ix]) - atan2(e1[ix],e1[ix]);
               }
	}
	gettimeofday(&endtime, NULL);
        te0=((double)(endtime.tv_sec*1000000-starttime.tv_sec*1000000+endtime.tv_usec-starttime.tv_usec))/1000000;
	te2 = 7 * ae * x100 / te0;
        printf("\n\tTest #2:\t Math.h trig:\t%lf/s", te2);
		/* Hyperbolic functions */
	gettimeofday(&starttime, NULL);
	for (i=0; i<x100; i++) {
#ifdef enable_threading
#pragma omp parallel for
#endif
		for (ix=0; ix<(ae-1); ix++)
		{
                e1[ix] = cosh(e1[ix]) + sinh(e1[ix]) - tanh(e1[ix]) * cosh(e1[ix]) / sinh(e1[ix]);
              	}	
	}
	gettimeofday(&endtime, NULL);
        te0=((double)(endtime.tv_sec*1000000-starttime.tv_sec*1000000+endtime.tv_usec-starttime.tv_usec))/1000000;
	te3 = 5 * ae * x100 / te0;
        printf("\n\tTest #3:\t Math.h hyperbolic:\t%lf/s", te3);
                /* Exponential and logarithmic functions */
	gettimeofday(&starttime, NULL);
	for (i=0; i<x100; i++) {
#ifdef enable_threading
#pragma omp parallel for
#endif
	        for (ix=0; ix<(ae-1); ix++)
                {
                e1[ix] = frexp(e1[ix], &n) + log(e1[ix]) - exp(e1[ix]) * ldexp(e1[ix], n) / log10(e1[ix]) + modf(e1[ix], &ip);
                }
	}
	gettimeofday(&endtime, NULL);
        te0=((double)(endtime.tv_sec*1000000-starttime.tv_sec*1000000+endtime.tv_usec-starttime.tv_usec))/1000000;
	te4 = 6 * ae * x100 / te0;
        printf("\n\tTest #4:\t Math.h exp/log:\t%lf/s", te4);
                /* Power functions */
	gettimeofday(&starttime, NULL);
	for (i=0; i<x100; i++) {
#ifdef enable_threading
#pragma omp parallel for
#endif
	        for (ix=0; ix<(ae-1); ix++)
                {
                e1[ix] = pow(e1[ix], e1[ix]) - sqrt(e1[ix]) - pow(e1[ix], e1[ix]) * sqrt(e1[ix]) / pow(e1[ix], e1[ix]) + sqrt(e1[ix]);
                }
	}
	gettimeofday(&endtime, NULL);
        te0=((double)(endtime.tv_sec*1000000-starttime.tv_sec*1000000+endtime.tv_usec-starttime.tv_usec))/1000000;
	te5 = 6 * ae * x100 / te0;
        printf("\n\tTest #5:\t Math.h power functions:\t%lf/s", te5);
                /* Rounding, absolute value and remainder functions */
	gettimeofday(&starttime, NULL);
	for (i=0; i<x100; i++) {
#ifdef enable_threading
#pragma omp parallel for
#endif
	        for (ix=0; ix<(ae-1); ix++)
                {
                e1[ix] = fmod(e1[ix], e1[ix]) - ceil(e1[ix]) * fabs(e1[ix]) / floor(e1[ix]);
                }
	}
	gettimeofday(&endtime, NULL);
	te0=((double)(endtime.tv_sec*1000000-starttime.tv_sec*1000000+endtime.tv_usec-starttime.tv_usec))/1000000;
	te6 = 4 * ae * x100 / te0;
	printf("\n\tTest #6:\t Math.h abs value, rounding, remainder:\t%lf/s\n", te6);
	pass0 =  (te1 + te2 + te3 + te4 + te5 + te6) / 6;
return pass0;
}
