parallel trapezoidal rule in C
0
I wrote a program to calculate the given integral as a function by extended trapezoidal rule. it seems that I go into an infinite loop but I can't understand why. I think the problem is in do_cmd function. or something with reference and pointers. please help. I try to implement the collected communication.
Thank you in advance.
#include<stdio.h>
#include <math.h>
#include "mpi.h"
int rank, size;
double part_of_integral(int n_start, int n_end, double a, double b, double h, double(*f)(double)){
double temp = f(n_start * h + a) + f(n_end * h + a);
temp *= h;
temp /= 2.0;
int i;
for (i = n_start + 1; i < n_end; ++i) {
temp += h * f(i * h);
}
return temp;
}
double func(double x){
return 1.0 / (1.0 + pow(x, 2.0));
}
double do_cmd(double* a, double* b, int* n){
MPI_Barrier(MPI_COMM_WORLD);
printf("a = %f, b = %f, n = %d, rank = %dn", *a, *b, *n, rank);
MPI_Bcast(a, 1, MPI_DOUBLE, 0, MPI_COMM_WORLD);
MPI_Bcast(b, 1, MPI_DOUBLE, 0, MPI_COMM_WORLD);
MPI_Bcast(n, 1, MPI_INT, 0, MPI_COMM_WORLD);
printf("aa= %f, bb = %f, nn = %d, rankb = %dn", *a, *b, *n, rank);
double h = (*b - *a) / *n;
int number_of_nodes = *n / size;
int n_start = rank * number_of_nodes;
int n_end = (rank + 1) * number_of_nodes;
double integral_output = part_of_integral(n_start, n_end, *a, *b, h, func);
double integral_result = 0.0;
MPI_Reduce(&integral_output, &integral_result, 1, MPI_DOUBLE, MPI_SUM, 0, MPI_COMM_WORLD);
if (rank == 0) {
integral_result *= 4.0;
}
return integral_result;
}
int main (int argc, char* argv) {
double a, b, tn, t2n;
int n;
double error = 1.0;
const double min_error = 1.0e-8;
MPI_Init(&argc, &argv);
MPI_Comm_size(MPI_COMM_WORLD, &size);
MPI_Comm_rank(MPI_COMM_WORLD, &rank);
tn = 0.0;
if (rank == 0) {
a = 0.0;
b = 1.0;
n = 2;
}
while(error > min_error) {
t2n = do_cmd(&a, &b, &n);
error = fabs(((3.0/4.0) * t2n) - ((1.0/3.0) * tn));
if (rank == 0) {
printf ("error: %fn", error);
n *= 2;
}
t2n = tn;
}
if (rank == 0) {
printf ("n = %d, 2n = %d, error = %fn", n/2, n, error);
}
MPI_Finalize();
return 0;
}
c function pointers parallel-processing
asked Nov 14 '18 at 21:43
DariyoushDariyoush
207
|
show 3 more comments
0
I wrote a program to calculate the given integral as a function by extended trapezoidal rule. it seems that I go into an infinite loop but I can't understand why. I think the problem is in do_cmd function. or something with reference and pointers. please help. I try to implement the collected communication.
Thank you in advance.
#include<stdio.h>
#include <math.h>
#include "mpi.h"
int rank, size;
double part_of_integral(int n_start, int n_end, double a, double b, double h, double(*f)(double)){
double temp = f(n_start * h + a) + f(n_end * h + a);
temp *= h;
temp /= 2.0;
int i;
for (i = n_start + 1; i < n_end; ++i) {
temp += h * f(i * h);
}
return temp;
}
double func(double x){
return 1.0 / (1.0 + pow(x, 2.0));
}
double do_cmd(double* a, double* b, int* n){
MPI_Barrier(MPI_COMM_WORLD);
printf("a = %f, b = %f, n = %d, rank = %dn", *a, *b, *n, rank);
MPI_Bcast(a, 1, MPI_DOUBLE, 0, MPI_COMM_WORLD);
MPI_Bcast(b, 1, MPI_DOUBLE, 0, MPI_COMM_WORLD);
MPI_Bcast(n, 1, MPI_INT, 0, MPI_COMM_WORLD);
printf("aa= %f, bb = %f, nn = %d, rankb = %dn", *a, *b, *n, rank);
double h = (*b - *a) / *n;
int number_of_nodes = *n / size;
int n_start = rank * number_of_nodes;
int n_end = (rank + 1) * number_of_nodes;
double integral_output = part_of_integral(n_start, n_end, *a, *b, h, func);
double integral_result = 0.0;
MPI_Reduce(&integral_output, &integral_result, 1, MPI_DOUBLE, MPI_SUM, 0, MPI_COMM_WORLD);
if (rank == 0) {
integral_result *= 4.0;
}
return integral_result;
}
int main (int argc, char* argv) {
double a, b, tn, t2n;
int n;
double error = 1.0;
const double min_error = 1.0e-8;
MPI_Init(&argc, &argv);
MPI_Comm_size(MPI_COMM_WORLD, &size);
MPI_Comm_rank(MPI_COMM_WORLD, &rank);
tn = 0.0;
if (rank == 0) {
a = 0.0;
b = 1.0;
n = 2;
}
while(error > min_error) {
t2n = do_cmd(&a, &b, &n);
error = fabs(((3.0/4.0) * t2n) - ((1.0/3.0) * tn));
if (rank == 0) {
printf ("error: %fn", error);
n *= 2;
}
t2n = tn;
}
if (rank == 0) {
printf ("n = %d, 2n = %d, error = %fn", n/2, n, error);
}
MPI_Finalize();
return 0;
}
c function pointers parallel-processing
asked Nov 14 '18 at 21:43
DariyoushDariyoush
207
1
Maybe the infinite loop is happening because of this:while(error > min_error).
– Fiddling Bits
Nov 14 '18 at 21:45
3
Did you mean to writetn = t2n? I don't seetngetting updated.
– David Cullen
Nov 14 '18 at 21:46
Addprintf("%f > %fn", error, min_error);belowerror = fabs(((3.0/4.0) * t2n) - ((1.0/3.0) * tn));.
– Fiddling Bits
Nov 14 '18 at 21:46
1
Inpart_of_integral, you are calling a functionf(which is not shown). Did you meanfunc? If so, the definition offunccomes too late. Compile with-Wall -Wextraand fix any warnings. I'd try to get this working with asizeof 1 first to be sure your calculation code and functions are correct. Then, increase the size. This separates the calc debug from debugging yourMPIcode.
– Craig Estey
Nov 14 '18 at 22:42
1
@CraigEsteyfis a pointer tofunc.
– David Cullen
Nov 14 '18 at 22:47
|
show 3 more comments
0
0
0
I wrote a program to calculate the given integral as a function by extended trapezoidal rule. it seems that I go into an infinite loop but I can't understand why. I think the problem is in do_cmd function. or something with reference and pointers. please help. I try to implement the collected communication.
Thank you in advance.
#include<stdio.h>
#include <math.h>
#include "mpi.h"
int rank, size;
double part_of_integral(int n_start, int n_end, double a, double b, double h, double(*f)(double)){
double temp = f(n_start * h + a) + f(n_end * h + a);
temp *= h;
temp /= 2.0;
int i;
for (i = n_start + 1; i < n_end; ++i) {
temp += h * f(i * h);
}
return temp;
}
double func(double x){
return 1.0 / (1.0 + pow(x, 2.0));
}
double do_cmd(double* a, double* b, int* n){
MPI_Barrier(MPI_COMM_WORLD);
printf("a = %f, b = %f, n = %d, rank = %dn", *a, *b, *n, rank);
MPI_Bcast(a, 1, MPI_DOUBLE, 0, MPI_COMM_WORLD);
MPI_Bcast(b, 1, MPI_DOUBLE, 0, MPI_COMM_WORLD);
MPI_Bcast(n, 1, MPI_INT, 0, MPI_COMM_WORLD);
printf("aa= %f, bb = %f, nn = %d, rankb = %dn", *a, *b, *n, rank);
double h = (*b - *a) / *n;
int number_of_nodes = *n / size;
int n_start = rank * number_of_nodes;
int n_end = (rank + 1) * number_of_nodes;
double integral_output = part_of_integral(n_start, n_end, *a, *b, h, func);
double integral_result = 0.0;
MPI_Reduce(&integral_output, &integral_result, 1, MPI_DOUBLE, MPI_SUM, 0, MPI_COMM_WORLD);
if (rank == 0) {
integral_result *= 4.0;
}
return integral_result;
}
int main (int argc, char* argv) {
double a, b, tn, t2n;
int n;
double error = 1.0;
const double min_error = 1.0e-8;
MPI_Init(&argc, &argv);
MPI_Comm_size(MPI_COMM_WORLD, &size);
MPI_Comm_rank(MPI_COMM_WORLD, &rank);
tn = 0.0;
if (rank == 0) {
a = 0.0;
b = 1.0;
n = 2;
}
while(error > min_error) {
t2n = do_cmd(&a, &b, &n);
error = fabs(((3.0/4.0) * t2n) - ((1.0/3.0) * tn));
if (rank == 0) {
printf ("error: %fn", error);
n *= 2;
}
t2n = tn;
}
if (rank == 0) {
printf ("n = %d, 2n = %d, error = %fn", n/2, n, error);
}
MPI_Finalize();
return 0;
}
c function pointers parallel-processing
asked Nov 14 '18 at 21:43
DariyoushDariyoush
207
I wrote a program to calculate the given integral as a function by extended trapezoidal rule. it seems that I go into an infinite loop but I can't understand why. I think the problem is in do_cmd function. or something with reference and pointers. please help. I try to implement the collected communication.
Thank you in advance.
#include<stdio.h>
#include <math.h>
#include "mpi.h"
int rank, size;
double part_of_integral(int n_start, int n_end, double a, double b, double h, double(*f)(double)){
double temp = f(n_start * h + a) + f(n_end * h + a);
temp *= h;
temp /= 2.0;
int i;
for (i = n_start + 1; i < n_end; ++i) {
temp += h * f(i * h);
}
return temp;
}
double func(double x){
return 1.0 / (1.0 + pow(x, 2.0));
}
double do_cmd(double* a, double* b, int* n){
MPI_Barrier(MPI_COMM_WORLD);
printf("a = %f, b = %f, n = %d, rank = %dn", *a, *b, *n, rank);
MPI_Bcast(a, 1, MPI_DOUBLE, 0, MPI_COMM_WORLD);
MPI_Bcast(b, 1, MPI_DOUBLE, 0, MPI_COMM_WORLD);
MPI_Bcast(n, 1, MPI_INT, 0, MPI_COMM_WORLD);
printf("aa= %f, bb = %f, nn = %d, rankb = %dn", *a, *b, *n, rank);
double h = (*b - *a) / *n;
int number_of_nodes = *n / size;
int n_start = rank * number_of_nodes;
int n_end = (rank + 1) * number_of_nodes;
double integral_output = part_of_integral(n_start, n_end, *a, *b, h, func);
double integral_result = 0.0;
MPI_Reduce(&integral_output, &integral_result, 1, MPI_DOUBLE, MPI_SUM, 0, MPI_COMM_WORLD);
if (rank == 0) {
integral_result *= 4.0;
}
return integral_result;
}
int main (int argc, char* argv) {
double a, b, tn, t2n;
int n;
double error = 1.0;
const double min_error = 1.0e-8;
MPI_Init(&argc, &argv);
MPI_Comm_size(MPI_COMM_WORLD, &size);
MPI_Comm_rank(MPI_COMM_WORLD, &rank);
tn = 0.0;
if (rank == 0) {
a = 0.0;
b = 1.0;
n = 2;
}
while(error > min_error) {
t2n = do_cmd(&a, &b, &n);
error = fabs(((3.0/4.0) * t2n) - ((1.0/3.0) * tn));
if (rank == 0) {
printf ("error: %fn", error);
n *= 2;
}
t2n = tn;
}
if (rank == 0) {
printf ("n = %d, 2n = %d, error = %fn", n/2, n, error);
}
MPI_Finalize();
return 0;
}
c function pointers parallel-processing
c function pointers parallel-processing
asked Nov 14 '18 at 21:43
DariyoushDariyoush
207
asked Nov 14 '18 at 21:43
DariyoushDariyoush
207
asked Nov 14 '18 at 21:43
DariyoushDariyoush
207
asked Nov 14 '18 at 21:43
DariyoushDariyoush
207
asked Nov 14 '18 at 21:43
DariyoushDariyoush
207
207
1
Maybe the infinite loop is happening because of this:while(error > min_error).
– Fiddling Bits
Nov 14 '18 at 21:45
3
Did you mean to writetn = t2n? I don't seetngetting updated.
– David Cullen
Nov 14 '18 at 21:46
Addprintf("%f > %fn", error, min_error);belowerror = fabs(((3.0/4.0) * t2n) - ((1.0/3.0) * tn));.
– Fiddling Bits
Nov 14 '18 at 21:46
1
Inpart_of_integral, you are calling a functionf(which is not shown). Did you meanfunc? If so, the definition offunccomes too late. Compile with-Wall -Wextraand fix any warnings. I'd try to get this working with asizeof 1 first to be sure your calculation code and functions are correct. Then, increase the size. This separates the calc debug from debugging yourMPIcode.
– Craig Estey
Nov 14 '18 at 22:42
1
@CraigEsteyfis a pointer tofunc.
– David Cullen
Nov 14 '18 at 22:47
|
show 3 more comments
1
Maybe the infinite loop is happening because of this:while(error > min_error).
– Fiddling Bits
Nov 14 '18 at 21:45
3
Did you mean to writetn = t2n? I don't seetngetting updated.
– David Cullen
Nov 14 '18 at 21:46
Addprintf("%f > %fn", error, min_error);belowerror = fabs(((3.0/4.0) * t2n) - ((1.0/3.0) * tn));.
– Fiddling Bits
Nov 14 '18 at 21:46
1
Inpart_of_integral, you are calling a functionf(which is not shown). Did you meanfunc? If so, the definition offunccomes too late. Compile with-Wall -Wextraand fix any warnings. I'd try to get this working with asizeof 1 first to be sure your calculation code and functions are correct. Then, increase the size. This separates the calc debug from debugging yourMPIcode.
– Craig Estey
Nov 14 '18 at 22:42
1
@CraigEsteyfis a pointer tofunc.
– David Cullen
Nov 14 '18 at 22:47
1
1
Maybe the infinite loop is happening because of this:
while(error > min_error).– Fiddling Bits
Nov 14 '18 at 21:45
Maybe the infinite loop is happening because of this:
while(error > min_error).– Fiddling Bits
Nov 14 '18 at 21:45
3
3
Did you mean to write
tn = t2n? I don't see tn getting updated.– David Cullen
Nov 14 '18 at 21:46
Did you mean to write
tn = t2n? I don't see tn getting updated.– David Cullen
Nov 14 '18 at 21:46
Add
printf("%f > %fn", error, min_error); below error = fabs(((3.0/4.0) * t2n) - ((1.0/3.0) * tn));.– Fiddling Bits
Nov 14 '18 at 21:46
Add
printf("%f > %fn", error, min_error); below error = fabs(((3.0/4.0) * t2n) - ((1.0/3.0) * tn));.– Fiddling Bits
Nov 14 '18 at 21:46
1
1
In
part_of_integral, you are calling a function f (which is not shown). Did you mean func? If so, the definition of func comes too late. Compile with -Wall -Wextra and fix any warnings. I'd try to get this working with a size of 1 first to be sure your calculation code and functions are correct. Then, increase the size. This separates the calc debug from debugging your MPI code.– Craig Estey
Nov 14 '18 at 22:42
In
part_of_integral, you are calling a function f (which is not shown). Did you mean func? If so, the definition of func comes too late. Compile with -Wall -Wextra and fix any warnings. I'd try to get this working with a size of 1 first to be sure your calculation code and functions are correct. Then, increase the size. This separates the calc debug from debugging your MPI code.– Craig Estey
Nov 14 '18 at 22:42
1
1
@CraigEstey
f is a pointer to func.– David Cullen
Nov 14 '18 at 22:47
@CraigEstey
f is a pointer to func.– David Cullen
Nov 14 '18 at 22:47
|
show 3 more comments
0
active
oldest
votes
Your Answer
StackExchange.ifUsing("editor", function () {
StackExchange.using("externalEditor", function () {
StackExchange.using("snippets", function () {
StackExchange.snippets.init();
});
});
}, "code-snippets");
StackExchange.ready(function() {
var channelOptions = {
tags: "".split(" "),
id: "1"
};
initTagRenderer("".split(" "), "".split(" "), channelOptions);
StackExchange.using("externalEditor", function() {
// Have to fire editor after snippets, if snippets enabled
if (StackExchange.settings.snippets.snippetsEnabled) {
StackExchange.using("snippets", function() {
createEditor();
});
}
else {
createEditor();
}
});
function createEditor() {
StackExchange.prepareEditor({
heartbeatType: 'answer',
autoActivateHeartbeat: false,
convertImagesToLinks: true,
noModals: true,
showLowRepImageUploadWarning: true,
reputationToPostImages: 10,
bindNavPrevention: true,
postfix: "",
imageUploader: {
brandingHtml: "Powered by u003ca class="icon-imgur-white" href="https://imgur.com/"u003eu003c/au003e",
contentPolicyHtml: "User contributions licensed under u003ca href="https://creativecommons.org/licenses/by-sa/3.0/"u003ecc by-sa 3.0 with attribution requiredu003c/au003e u003ca href="https://stackoverflow.com/legal/content-policy"u003e(content policy)u003c/au003e",
allowUrls: true
},
onDemand: true,
discardSelector: ".discard-answer"
,immediatelyShowMarkdownHelp:true
});
}
});
draft saved
draft discarded
Sign up or log in
StackExchange.ready(function () {
StackExchange.helpers.onClickDraftSave('#login-link');
});
Sign up using Google
Sign up using Facebook
Sign up using Email and Password
Post as a guest
Required, but never shown
StackExchange.ready(
function () {
StackExchange.openid.initPostLogin('.new-post-login', 'https%3a%2f%2fstackoverflow.com%2fquestions%2f53309169%2fparallel-trapezoidal-rule-in-c%23new-answer', 'question_page');
}
);
Post as a guest
Required, but never shown
0
active
oldest
votes
0
active
oldest
votes
active
oldest
votes
active
oldest
votes
draft saved
draft discarded
Thanks for contributing an answer to Stack Overflow!
- Please be sure to answer the question. Provide details and share your research!
But avoid …
- Asking for help, clarification, or responding to other answers.
- Making statements based on opinion; back them up with references or personal experience.
To learn more, see our tips on writing great answers.
draft saved
draft discarded
Sign up or log in
StackExchange.ready(function () {
StackExchange.helpers.onClickDraftSave('#login-link');
});
Sign up using Google
Sign up using Facebook
Sign up using Email and Password
Post as a guest
Required, but never shown
StackExchange.ready(
function () {
StackExchange.openid.initPostLogin('.new-post-login', 'https%3a%2f%2fstackoverflow.com%2fquestions%2f53309169%2fparallel-trapezoidal-rule-in-c%23new-answer', 'question_page');
}
);
Post as a guest
Required, but never shown
Sign up or log in
StackExchange.ready(function () {
StackExchange.helpers.onClickDraftSave('#login-link');
});
Sign up using Google
Sign up using Facebook
Sign up using Email and Password
Post as a guest
Required, but never shown
Sign up or log in
StackExchange.ready(function () {
StackExchange.helpers.onClickDraftSave('#login-link');
});
Sign up using Google
Sign up using Facebook
Sign up using Email and Password
Post as a guest
Required, but never shown
Sign up or log in
StackExchange.ready(function () {
StackExchange.helpers.onClickDraftSave('#login-link');
});
Sign up using Google
Sign up using Facebook
Sign up using Email and Password
Sign up using Google
Sign up using Facebook
Sign up using Email and Password
Post as a guest
Required, but never shown
Required, but never shown
Required, but never shown
Required, but never shown
Required, but never shown
Required, but never shown
Required, but never shown
Required, but never shown
Required, but never shown
1
Maybe the infinite loop is happening because of this:
while(error > min_error).– Fiddling Bits
Nov 14 '18 at 21:45
3
Did you mean to write
tn = t2n? I don't seetngetting updated.– David Cullen
Nov 14 '18 at 21:46
Add
printf("%f > %fn", error, min_error);belowerror = fabs(((3.0/4.0) * t2n) - ((1.0/3.0) * tn));.– Fiddling Bits
Nov 14 '18 at 21:46
1
In
part_of_integral, you are calling a functionf(which is not shown). Did you meanfunc? If so, the definition offunccomes too late. Compile with-Wall -Wextraand fix any warnings. I'd try to get this working with asizeof 1 first to be sure your calculation code and functions are correct. Then, increase the size. This separates the calc debug from debugging yourMPIcode.– Craig Estey
Nov 14 '18 at 22:42
1
@CraigEstey
fis a pointer tofunc.– David Cullen
Nov 14 '18 at 22:47