RInside Version 0.2.12
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros
rinside_mpi_sample3.cpp
Go to the documentation of this file.
1 // -*- mode: C++; c-indent-level: 4; c-basic-offset: 4; tab-width: 8; -*-
2 //
3 // Simple mpi example: simulate sampling/averaging on multiple nodes and gathering the results.
4 //
5 // MPI C++ API version of file contributed by Jianping Hua
6 //
7 // Copyright (C) 2010 - 2011 Dirk Eddelbuettel and Romain Francois
8 //
9 // GPL'ed
10 
11 #include <mpi.h> // mpi header file
12 #include <RInside.h> // for the embedded R via RInside
13 
14 int main(int argc, char *argv[]) {
15 
16  MPI::Init(argc, argv); // mpi initialization
17  int myrank = MPI::COMM_WORLD.Get_rank(); // obtain current node rank
18  int nodesize = MPI::COMM_WORLD.Get_size(); // obtain total nodes running.
19 
20  int sndcnt = 1, rcvcnt = 1; // # of elements in send/recv buffer
21  double sendValue; // value to be collected in current node
22  double *allvalues = new double[nodesize]; // to save all results
23 
24  // simulation info
25  // to sample from a uniform distribution
26  int rangeMin = 0, rangeMax = 10; // range of uniform distribution
27  int sampleSize = 2; // points in each sample
28 
29  try {
30  RInside R(argc, argv); // create an embedded R instance
31 
32  std::stringstream txt;
33  txt << "x <- " << rangeMin << std::endl;
34  R.parseEvalQ( txt.str() ); // assign x with lower range of uniform distribution
35 
36  txt << "y <- " << rangeMax << std::endl;
37  R.parseEvalQ( txt.str() ); // assign y with upper range of uniform distribution
38 
39  txt << "n <- " << sampleSize << std::endl;
40  R.parseEvalQ( txt.str() ); // assign n with the size of sample
41 
42  std::string evalstr = "mean(runif(n,x,y))"; // sampling, compute the mean
43  Rcpp::NumericVector m = R.parseEval(evalstr); // eval str, convert result to NumericVector
44  sendValue = m( 0 ); // assign the return value to the variable to be gathered
45 
46  //gather together values from all processes to allvalues
47  MPI::COMM_WORLD.Gather((const void*)&sendValue, sndcnt, MPI::DOUBLE, (void*)allvalues, rcvcnt, MPI::DOUBLE, 0);
48 
49  // show what inidividual node's contribution
50  std::cout << "node " << myrank << " has mean " << m(0) << std::endl;
51 
52  } catch(std::exception& ex) {
53  std::cerr << "Exception caught: " << ex.what() << std::endl;
54  } catch(...) {
55  std::cerr << "Unknown exception caught" << std::endl;
56  }
57 
58  // show gathered results in node 0
59  if ( myrank == 0 ) {
60  std::cout << "values of all " << nodesize << " trials: " << std::endl;
61  for ( int i = 0; i < nodesize; i++ )
62  std::cout << allvalues[ i ] << ", ";
63  std::cout << std::endl;
64  }
65 
66  // clean up
67  delete[] allvalues;
68 
69  MPI::Finalize(); // mpi finalization
70 
71  exit(0);
72 }
73 
void parseEvalQ(const std::string &line)
Definition: RInside.cpp:366
int parseEval(const std::string &line, SEXP &ans)
Definition: RInside.cpp:308
int main(int argc, char *argv[])