-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy pathmpi_trap4.c
189 lines (161 loc) · 6 KB
/
mpi_trap4.c
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
/* File: mpi_trap4.c
* Purpose: Use MPI to implement a parallel version of the trapezoidal
* rule. This version uses collective communications and
* MPI derived datatypes to distribute the input data and
* compute the global sum.
*
* Input: The endpoints of the interval of integration and the number
* of trapezoids
* Output: Estimate of the integral from a to b of f(x)
* using the trapezoidal rule and n trapezoids.
*
* Compile: mpicc -g -Wall -o mpi_trap4 mpi_trap4.c
* Run: mpiexec -n <number of processes> ./mpi_trap4
*
* Algorithm:
* 1. Each process calculates "its" interval of
* integration.
* 2. Each process estimates the integral of f(x)
* over its interval using the trapezoidal rule.
* 3a. Each process != 0 sends its integral to 0.
* 3b. Process 0 sums the calculations received from
* the individual processes and prints the result.
*
* Note: f(x) is all hardwired.
*
* IPP: Section 3.5 (pp. 117 and ff.)
*/
#include <stdio.h>
/* We'll be using MPI routines, definitions, etc. */
#include <mpi.h>
/* Build a derived datatype for distributing the input data */
void Build_mpi_type(double* a_p, double* b_p, int* n_p,
MPI_Datatype* input_mpi_t_p);
/* Get the input values */
void Get_input(int my_rank, int comm_sz, double* a_p, double* b_p,
int* n_p);
/* Calculate local integral */
double Trap(double left_endpt, double right_endpt, int trap_count,
double base_len);
/* Function we're integrating */
double f(double x);
int main(void) {
int my_rank, comm_sz, n, local_n;
double a, b, h, local_a, local_b;
double local_int, total_int;
/* Let the system do what it needs to start up MPI */
MPI_Init(NULL, NULL);
/* Get my process rank */
MPI_Comm_rank(MPI_COMM_WORLD, &my_rank);
/* Find out how many processes are being used */
MPI_Comm_size(MPI_COMM_WORLD, &comm_sz);
Get_input(my_rank, comm_sz, &a, &b, &n);
h = (b-a)/n; /* h is the same for all processes */
local_n = n/comm_sz; /* So is the number of trapezoids */
/* Length of each process' interval of
* integration = local_n*h. So my interval
* starts at: */
local_a = a + my_rank*local_n*h;
local_b = local_a + local_n*h;
local_int = Trap(local_a, local_b, local_n, h);
/* Add up the integrals calculated by each process */
MPI_Reduce(&local_int, &total_int, 1, MPI_DOUBLE, MPI_SUM, 0,
MPI_COMM_WORLD);
/* Print the result */
if (my_rank == 0) {
printf("With n = %d trapezoids, our estimate\n", n);
printf("of the integral from %f to %f = %.15e\n",
a, b, total_int);
}
/* Shut down MPI */
MPI_Finalize();
return 0;
} /* main */
/*------------------------------------------------------------------
* Function: Build_mpi_type
* Purpose: Build a derived datatype so that the three
* input values can be sent in a single message.
* Input args: a_p: pointer to left endpoint
* b_p: pointer to right endpoint
* n_p: pointer to number of trapezoids
* Output args: input_mpi_t_p: the new MPI datatype
*/
void Build_mpi_type(
double* a_p /* in */,
double* b_p /* in */,
int* n_p /* in */,
MPI_Datatype* input_mpi_t_p /* out */) {
int array_of_blocklengths[3] = {1, 1, 1};
MPI_Datatype array_of_types[3] = {MPI_DOUBLE, MPI_DOUBLE, MPI_INT};
MPI_Aint a_addr, b_addr, n_addr;
MPI_Aint array_of_displacements[3] = {0};
MPI_Get_address(a_p, &a_addr);
MPI_Get_address(b_p, &b_addr);
MPI_Get_address(n_p, &n_addr);
array_of_displacements[1] = b_addr-a_addr;
array_of_displacements[2] = n_addr-a_addr;
MPI_Type_create_struct(3, array_of_blocklengths,
array_of_displacements, array_of_types,
input_mpi_t_p);
MPI_Type_commit(input_mpi_t_p);
} /* Build_mpi_type */
/*------------------------------------------------------------------
* Function: Get_input
* Purpose: Get the user input: the left and right endpoints
* and the number of trapezoids
* Input args: my_rank: process rank in MPI_COMM_WORLD
* comm_sz: number of processes in MPI_COMM_WORLD
* Output args: a_p: pointer to left endpoint
* b_p: pointer to right endpoint
* n_p: pointer to number of trapezoids
*/
void Get_input(
int my_rank /* in */,
int comm_sz /* in */,
double* a_p /* out */,
double* b_p /* out */,
int* n_p /* out */) {
MPI_Datatype input_mpi_t;
Build_mpi_type(a_p, b_p, n_p, &input_mpi_t);
if (my_rank == 0) {
printf("Enter a, b, and n\n");
scanf("%lf %lf %d", a_p, b_p, n_p);
}
MPI_Bcast(a_p, 1, input_mpi_t, 0, MPI_COMM_WORLD);
MPI_Type_free(&input_mpi_t);
} /* Get_input */
/*------------------------------------------------------------------
* Function: Trap
* Purpose: Serial function for estimating a definite integral
* using the trapezoidal rule
* Input args: left_endpt
* right_endpt
* trap_count
* base_len
* Return val: Trapezoidal rule estimate of integral from
* left_endpt to right_endpt using trap_count
* trapezoids
*/
double Trap(
double left_endpt /* in */,
double right_endpt /* in */,
int trap_count /* in */,
double base_len /* in */) {
double estimate, x;
int i;
estimate = (f(left_endpt) + f(right_endpt))/2.0;
for (i = 1; i <= trap_count-1; i++) {
x = left_endpt + i*base_len;
estimate += f(x);
}
estimate = estimate*base_len;
return estimate;
} /* Trap */
/*------------------------------------------------------------------
* Function: f
* Purpose: Compute value of function to be integrated
* Input args: x
*/
double f(double x /* in */) {
return x*x;
} /* f */