RInside Version 0.2.12
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros
rinside_mpi_sample1.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 // This file was contributed by Jianping Hua
6 //
7 // Copyright (C) 2010 - 2011 Jianping Hua, 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 related
17  int myrank, nodesize; // node information
18  int sndcnt = 1, rcvcnt = 1; // # of elements in send/recv buffer
19  MPI_Init(&argc,&argv); // mpi initialization
20  MPI_Comm_rank(MPI_COMM_WORLD, &myrank); // obtain current node rank
21  MPI_Comm_size(MPI_COMM_WORLD, &nodesize); // obtain total nodes running.
22  double sendValue; // value to be collected in current node
23  double *allvalues = new double[nodesize]; // to save all results
24 
25  // simulation info
26  // to sample from a uniform distribution
27  int rangeMin = 0, rangeMax = 10; // range of uniform distribution
28  int sampleSize = 2; // points in each sample
29 
30  try {
31  RInside R(argc, argv); // create an embedded R instance
32 
33  std::stringstream txt;
34  txt << "x <- " << rangeMin << std::endl;
35  R.parseEvalQ( txt.str() ); // assign x with lower range of uniform distribution
36 
37  txt << "y <- " << rangeMax << std::endl;
38  R.parseEvalQ( txt.str() ); // assign y with upper range of uniform distribution
39 
40  txt << "n <- " << sampleSize << std::endl;
41  R.parseEvalQ( txt.str() ); // assign n with the size of sample
42 
43  std::string evalstr = "mean(runif(n,x,y))"; // sampling, compute the mean
44  Rcpp::NumericVector m = R.parseEval(evalstr); // eval str, convert result to NumericVector
45  sendValue = m( 0 ); // assign the return value to the variable to be gathered
46 
47  //gather together values from all processes to allvalues
48  MPI_Gather(&sendValue, sndcnt, MPI_DOUBLE, allvalues, rcvcnt, MPI_DOUBLE, 0, MPI_COMM_WORLD);
49 
50  // show what inidividual node's contribution
51  std::cout << "node " << myrank << " has mean " << m(0) << std::endl;
52 
53  } catch(std::exception& ex) {
54  std::cerr << "Exception caught: " << ex.what() << std::endl;
55  } catch(...) {
56  std::cerr << "Unknown exception caught" << std::endl;
57  }
58 
59  // show gathered results in node 0
60  if ( myrank == 0 ) {
61  std::cout << "values of all " << nodesize << " trials: " << std::endl;
62  for ( int i = 0; i < nodesize; i++ )
63  std::cout << allvalues[ i ] << ", ";
64  std::cout << std::endl;
65  }
66 
67  // clean up
68  delete[] allvalues;
69 
70  MPI_Finalize(); // mpi finalization
71 
72  exit(0);
73 }
74 
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[])