ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/OpenMD/trunk/test/math/testRandNumGen.cpp
(Generate patch)

Comparing trunk/test/math/testRandNumGen.cpp (file contents):
Revision 391 by tim, Wed Mar 2 07:32:15 2005 UTC vs.
Revision 1390 by gezelter, Wed Nov 25 20:02:06 2009 UTC

# Line 1 | Line 1
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 OpenMD;
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) {
# Line 70 | Line 81 | void RandNumGenTestCase::testParallelRandNumGen(){    
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  
# Line 78 | Line 92 | void RandNumGenTestCase::testParallelRandNumGen(){    
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          }
# Line 86 | Line 100 | void RandNumGenTestCase::testParallelRandNumGen(){    
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();
# Line 100 | Line 113 | int main(int argc, char* argv[]) {
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   }

Diff Legend

– Removed lines
+ Added lines
< Changed lines
> Changed lines