Frage

I am trying to generate some uniform real numbers for a Monte Carlo integration but the routine I build was returning some really strange values. Upon closer inspection I notices that Boost was returning some crazy looking random numbers e.g.:

     temp = -0.185276
     temp = -0.864523
     temp = -0.0942081
     temp = -0.164991
     temp = -0.873013
     temp = -0.0311322
     temp = -0.0866241
     temp = -0.778966
     temp = -0.367641
     temp = -0.691833
     temp = 5.66499e-310
     temp = 9.42007e-311
     temp = 6.29821e-310
     temp = 5.80603e-310
     temp = 8.82973e-311
     temp = 6.73679e-310
     temp = 6.35094e-310
     temp = 1.53691e-310
     temp = 4.39696e-310
     temp = 2.14277e-310

Whilst these numbers are technically still reals generated between the bounds -1 and 1 I would prefer it if they weren't quite so small!

My implementation of the call to boost is in a function which is called multiple times (for different bounding values) as follows:

// Define Boost typedefs
typedef boost::mt19937 Engine;
typedef boost::uniform_real<double> Distribution;
typedef boost::variate_generator <Engine, Distribution> Generator;

int main (void) {
    ...
    Integral = MCRecursion(...);
    ...
    return 0;
}

double MCRecursion (int Count, double Lower, double Upper, double (*Integrand)(double)) {

// Define Boost objects
Engine       Eng;
Distribution Dist (Lower, Upper);
Generator    RandomGen (Eng, Dist);

Eng.seed(time(0));

// Variables for Monte Carlo sample sums
double Sum = 0.0;
double temp;

for (int i = 0; i < Count; i++) {
    temp = RandomGen();
    std::cout << "         temp = " << temp << std::endl;
    Sum += Integrand(temp);
}

return (Upper - Lower) * Sum / Count;
}

I assume the problem is something with my implementation but I can't find any errors. Any and all help appreciated! Cheers, Jack

EDIT

Code for calling MCRecursion:

The Code I am writting runs a Monte Carlo on the entire domain I am interested in [Lower, Upper] and then looks again at the left half of the whole domain and the right half of the domain. e.g. if we were integrating f(x) between -a and a I calculate the full integral using:

double FullDomain = MCRecursion (1e5, LowerBound, UpperBound, f);
double Centre = (Upper + Lower) / 2.0;
double LeftHalf = MCRecursion (1e5, LowerBound, Centre, f);
double RightHalf = MCRecursion (1e5, Centre, UpperBound, f);

and I then look at the uncertainty by calculating: double difference = fabs(FullDomain - LeftHalf - Righthalf); to see if more samples is 'worth it' in some sense Jack

War es hilfreich?

Lösung

Based on the pastebin the questioner posted in the comments:

This is not a problem with the random library but rather a simple programming error. Compiling the code throws the warning:

../src/Test.cpp: In function ‘double AdaptiveMCRecursion(double, double, double, double, int, double, double (*)(double))’:
../src/Test.cpp:100:72: warning: ‘Right’ is used uninitialized in this function [-Wuninitialized]
         double Right = MCSample (Count, Central, Right,   Integrand);

So all the behaviour from that line on is basically undefined. Especially it results in calling the function MCSample with an undetermined Upper parameter. So your result is not unexpected. You are actually lucky the program runs at all.

Lizenziert unter: CC-BY-SA mit Zuschreibung
Nicht verbunden mit StackOverflow
scroll top