Cyclops Tensor Framework
parallel arithmetic on multidimensional arrays
endomorphism.cxx
Go to the documentation of this file.
1 /*Copyright (c) 2011, Edgar Solomonik, all rights reserved.*/
2 
10 #include <ctf.hpp>
11 using namespace CTF;
12 
13 void fdbl(double & a){
14  a=a*a*a;
15 }
16 
17 int endomorphism(int n,
18  World & dw){
19 
20  int shapeN4[] = {NS,NS,NS,NS};
21  int sizeN4[] = {n+1,n,n+2,n+3};
22 
23  Tensor<> A(4, sizeN4, shapeN4, dw);
24 
25  A.fill_random(-.5, .5);
26 
27 
28  double * all_start_data;
29  int64_t nall;
30  A.read_all(&nall, &all_start_data);
31 
32  double scale = 1.0;
33 
34  CTF::Transform<double> endo([=](double & d){ d=scale*d*d*d; });
35  // below is equivalent to A.scale(1.0, "ijkl", endo);
36  endo(A["ijkl"]);
37 
38  double * all_end_data;
39  int64_t nall2;
40  A.read_all(&nall2, &all_end_data);
41 
42  int pass = (nall == nall2);
43  if (pass){
44  for (int64_t i=0; i<nall; i++){
45  fdbl(all_start_data[i]);
46  if (fabs(all_start_data[i]-all_end_data[i])>=1.E-6) pass =0;
47  }
48  }
49  MPI_Allreduce(MPI_IN_PLACE, &pass, 1, MPI_INT, MPI_MIN, MPI_COMM_WORLD);
50 
51  if (dw.rank == 0){
52  if (pass){
53  printf("{ A[\"ijkl\"] = A[\"ijkl\"]^3 } passed\n");
54  } else {
55  printf("{ A[\"ijkl\"] = A[\"ijkl\"]^3 } failed\n");
56  }
57  }
58 
59  delete [] all_start_data;
60  delete [] all_end_data;
61 
62  return pass;
63 }
64 
65 
66 #ifndef TEST_SUITE
67 
68 char* getCmdOption(char ** begin,
69  char ** end,
70  const std::string & option){
71  char ** itr = std::find(begin, end, option);
72  if (itr != end && ++itr != end){
73  return *itr;
74  }
75  return 0;
76 }
77 
78 
79 int main(int argc, char ** argv){
80  int rank, np, n;
81  int const in_num = argc;
82  char ** input_str = argv;
83 
84  MPI_Init(&argc, &argv);
85  MPI_Comm_rank(MPI_COMM_WORLD, &rank);
86  MPI_Comm_size(MPI_COMM_WORLD, &np);
87 
88  if (getCmdOption(input_str, input_str+in_num, "-n")){
89  n = atoi(getCmdOption(input_str, input_str+in_num, "-n"));
90  if (n < 0) n = 5;
91  } else n = 5;
92 
93 
94  {
95  World dw(MPI_COMM_WORLD, argc, argv);
96 
97  if (rank == 0){
98  printf("Computing endomorphism A_ijkl = f(A_ijkl)\n");
99  }
100  endomorphism(n, dw);
101  }
102 
103 
104  MPI_Finalize();
105  return 0;
106 }
107 
113 #endif
int main(int argc, char **argv)
def rank(self)
Definition: core.pyx:312
void read_all(int64_t *npair, dtype **data, bool unpack=false)
collects the entire tensor data on each process (not memory scalable)
Definition: tensor.cxx:377
Definition: common.h:37
an instance of the CTF library (world) on a MPI communicator
Definition: world.h:19
string
Definition: core.pyx:456
char * getCmdOption(char **begin, char **end, const std::string &option)
int endomorphism(int n, World &dw)
void fill_random(dtype rmin, dtype rmax)
fills local unique tensor elements to random values in the range [min,max] works only for dtype in {f...
Definition: tensor.cxx:928
def scale(self, scl)
Definition: core.pyx:325
int rank
rank of local processor
Definition: world.h:24
Definition: apsp.cxx:17
an instance of a tensor within a CTF world
Definition: tensor.h:74
void fdbl(double &a)
def np(self)
Definition: core.pyx:315