/* File: omp_trap1a.c * Purpose: Estimate definite integral using trapezoidal * rule. * * Input: a, b, n * Output: estimate of integral from a to b of f(x) * using n trapezoids. * * Compile: gcc -g -Wall -fopenmp -o omp_trap1a omp_trap1a.c * Usage: ./omp_trap1a * * Notes: * 1. The function f(x) is hardwired. * 2. In this version, each thread explicitly computes the integral * over its assigned subinterval, a critical directive is used * for the global sum. * 3. This version assumes that n is evenly divisible by the * number of threads */ #include #include #include #include void Usage(char* prog_name); double f(double x); /* Function we're integrating */ void Trap(double a, double b, int n, double* integral_p); int main(int argc, char* argv[]) { double integral; /* Store result in integral */ double a, b; /* Left and right endpoints */ int n; /* Total number of trapezoids */ int thread_count; if (argc != 2) Usage(argv[0]); thread_count = strtol(argv[1], NULL, 10); printf("Enter a, b, and n\n"); scanf("%lf %lf %d", &a, &b, &n); integral = 0.0; # pragma omp parallel num_threads(thread_count) Trap(a, b, n, &integral); printf("With n = %d trapezoids, our estimate\n", n); printf("of the integral from %f to %f = %.15e\n", a, b, integral); return 0; } /* main */ /*-------------------------------------------------------------------- * Function: Usage * Purpose: Print command line for function and terminate * In arg: prog_name */ void Usage(char* prog_name) { fprintf(stderr, "usage: %s \n", prog_name); exit(0); } /* Usage */ /*------------------------------------------------------------------ * Function: f * Purpose: Compute value of function to be integrated * Input arg: x * Return val: f(x) */ double f(double x) { double return_val; return_val = x*x; return return_val; } /* f */ /*------------------------------------------------------------------ * Function: Trap * Purpose: Use trapezoidal rule to compute definite integral * Input args: * a: left endpoint * b: right endpoint * n: number of trapezoids * Return value: Estimate of Integral from a to b of f(x) */ void Trap(double a, double b, int n, double* integral_p) { double h, x, my_integral; double local_a, local_b; int i, local_n; int my_rank = omp_get_thread_num(); int thread_count = omp_get_num_threads(); h = (b-a)/n; local_n = n/thread_count; local_a = a + my_rank*local_n*h; local_b = local_a + local_n*h; my_integral = (f(local_a) + f(local_b))/2.0; for (i = 1; i <= local_n-1; i++) { x = local_a + i*h; my_integral += f(x); } my_integral = my_integral*h; # pragma omp critical *integral_p += my_integral; } /* Trap */