| 1 | + | #include <functional> | 
| 2 |  | #include <cassert> | 
| 3 |  | #include <fstream> | 
| 4 |  | #include <algorithm> | 
| 5 | < | #include "utils/simError.h" | 
| 5 | > | #include <vector> | 
| 6 | > | #include <iostream> | 
| 7 | > |  | 
| 8 | > | #include "math/SeqRandNumGen.hpp" | 
| 9 | > | #ifdef IS_MPI | 
| 10 | > | #include <mpi.h> | 
| 11 |  | #include "math/ParallelRandNumGen.hpp" | 
| 12 | + | #endif | 
| 13 | + | using namespace oopse; | 
| 14 | + | using namespace std; | 
| 15 |  |  | 
| 16 | < | void RandNumGenTestCase::testUniform(){ | 
| 17 | < | MTRand randNumGen(823645754); | 
| 18 | < |  | 
| 10 | < | const int N = 16; | 
| 16 | > | void testUniform(){ | 
| 17 | > | SeqRandNumGen randNumGen(823645754); | 
| 18 | > | const int N = 100; | 
| 19 |  | std::vector<unsigned long int> histogram(N, 0); | 
| 20 | < | const int num = 1000000; | 
| 20 | > | const unsigned long int num = 10000000; | 
| 21 |  | for (int i = 0; i <num; ++i) { | 
| 22 | < | ++histogram[randNumGen.randInt(N -1 )]; // rantInt returns an integer in [0, N-1] | 
| 22 | > | int index = randNumGen.randInt(N -1 ); | 
| 23 | > | ++histogram[index]; // rantInt returns an integer in [0, N-1] | 
| 24 |  | } | 
| 25 | < |  | 
| 26 | < | std::ofstream uniform("uniform.dat") | 
| 18 | < | int avg = num / N; | 
| 25 | > | std::ofstream uniform("uniform.dat"); | 
| 26 | > | double avg = num / N; | 
| 27 |  | double tolerance = 0.01*avg; | 
| 28 | < | for (int i = 0; i < num; ++i) { | 
| 29 | < | assert((histogram[i] - avg) /avg <= tolerance); | 
| 28 | > | for (int i = 0; i < N; ++i) { | 
| 29 | > | //assert((histogram[i] - avg) /avg <= tolerance); | 
| 30 |  | uniform << i << "\t" << histogram[i] << std::endl; | 
| 31 |  | } | 
| 32 |  | } | 
| 33 |  |  | 
| 34 | < | void RandNumGenTestCase::testGaussian(){ | 
| 35 | < | MTRand randNumGen(823645754); | 
| 34 | > | void testGaussian(){ | 
| 35 | > | SeqRandNumGen randNumGen(823645754); | 
| 36 |  | double mean = 100.0; | 
| 37 |  | double variance = 1.0; | 
| 38 | < | const int num = 1000000; | 
| 38 | > | const unsigned long int num = 1000000; | 
| 39 |  | double interval = 0.1; | 
| 40 |  | const int size = 2000; | 
| 41 |  | vector<unsigned long int> histogram(size , 0); | 
| 42 |  | vector<double> normalizedHistogram(size); | 
| 43 | < | for (int i = 0; i < num; ++i) { | 
| 43 | > | for (unsigned long int i = 0; i < num; ++i) { | 
| 44 |  | int index = static_cast<int>(randNumGen.randNorm(mean, variance) / interval); | 
| 45 |  | ++histogram[index]; | 
| 46 |  | } | 
| 47 |  |  | 
| 48 |  | std::transform(histogram.begin(), histogram.end(), normalizedHistogram.begin(), std::bind2nd(std::divides<double>(), num)); | 
| 49 |  | std::ofstream gaussian("gaussian.dat"); | 
| 50 | < | for (int i = 0; i < num; ++i) { | 
| 51 | < | gaussian << i << "\t" << normalizedHistogram[i] << std::endl; | 
| 50 | > | for (int i = 0; i < size; ++i) { | 
| 51 | > | gaussian << i*interval << "\t" << normalizedHistogram[i] << std::endl; | 
| 52 |  | } | 
| 53 |  | } | 
| 54 | < |  | 
| 55 | < | void RandNumGenTestCase::testParallelRandNumGen(){ | 
| 56 | < | const int seed = 324271632; | 
| 57 | < | const int nloops = 1000000; | 
| 54 | > | #ifdef IS_MPI | 
| 55 | > | void testParallelRandNumGen(){ | 
| 56 | > | const unsigned long int seed = 324271632; | 
| 57 | > | const unsigned long int nloops = 1000000; | 
| 58 |  | MPI_Status istatus; | 
| 59 |  | ParallelRandNumGen mpiRandNumGen(seed); | 
| 60 |  | const int masterNode = 0; | 
| 61 | < | if (worldRank = masterNode) { | 
| 61 | > | int myRank; | 
| 62 | > | MPI_Comm_rank( MPI_COMM_WORLD, &myRank ); | 
| 63 | > | if (myRank == masterNode) { | 
| 64 |  |  | 
| 65 | < | MTRand singleRandNumGen(seed); | 
| 66 | < |  | 
| 65 | > | SeqRandNumGen singleRandNumGen(seed); | 
| 66 | > | std::ofstream singleOs("single.dat"); | 
| 67 | > | std::ofstream parallelOs("parallel.dat"); | 
| 68 |  | int nProcessors; | 
| 69 |  | MPI_Comm_size(MPI_COMM_WORLD, &nProcessors); | 
| 70 |  | std::vector<unsigned long int> mpiRandNums(nProcessors); | 
| 71 |  | std::vector<unsigned long int> singleRandNums(nProcessors); | 
| 72 |  |  | 
| 73 | < | for (int i = 0; i < nloops; ++i) { | 
| 73 | > | for (unsigned long int i = 0; i < nloops; ++i) { | 
| 74 |  | mpiRandNums[masterNode] = mpiRandNumGen.randInt(); | 
| 75 |  |  | 
| 76 |  | for (int j = 0; j < nProcessors; ++j) { | 
| 81 |  | singleRandNums[j] = mpiRandNumGen.randInt(); | 
| 82 |  | } | 
| 83 |  |  | 
| 84 | < | assert(mpiRandNums, singleRandNums); | 
| 84 | > | for (int j = 0; j < nProcessors; ++j) { | 
| 85 | > | singleOs << singleRandNums[j] << "\n"; | 
| 86 | > | parallelOs << singleRandNums[j] << "\n"; | 
| 87 | > | } | 
| 88 |  | } | 
| 89 |  |  | 
| 90 |  |  | 
| 92 |  | } else { | 
| 93 |  |  | 
| 94 |  | unsigned long int randNum; | 
| 95 | < | for (int i = 0; i < nloops; ++i) { | 
| 95 | > | for (unsigned long int i = 0; i < nloops; ++i) { | 
| 96 |  | randNum = mpiRandNumGen.randInt(); | 
| 97 |  | MPI_Send(&randNum, 1, MPI_INT, masterNode, i, MPI_COMM_WORLD); | 
| 98 |  | } | 
| 100 |  | } | 
| 101 |  |  | 
| 102 |  | } | 
| 103 | < |  | 
| 90 | < |  | 
| 103 | > | #endif | 
| 104 |  | int main(int argc, char* argv[]) { | 
| 105 | < |  | 
| 106 | < | MPI_Init(argc, argv); | 
| 107 | < |  | 
| 105 | > | #ifdef IS_MPI | 
| 106 | > | MPI_Init(&argc, &argv); | 
| 107 | > | std::cout << "begin test" << std::endl; | 
| 108 |  | if (worldRank == 0 ) { | 
| 109 |  | testUniform(); | 
| 110 |  | testGaussian(); | 
| 113 |  | testParallelRandNumGen(); | 
| 114 |  |  | 
| 115 |  | MPI_Finalize(); | 
| 116 | + | #else | 
| 117 | + | std::cout << "begin test" <<std::endl; | 
| 118 | + | testUniform(); | 
| 119 | + | testGaussian(); | 
| 120 | + | #endif | 
| 121 | + | std::cout << "test done" << std::endl; | 
| 122 | + | return 0; | 
| 123 |  | } |