noob problem with a number generator

hey guys thanks for reading this.. i know that its a simple error on my part as i dont really know what im doing..

Am using this "Numerical recipes in C" to help me with some complex maths, but my programming skills are shot. The book can be found in web-form at http://www.ma.utexas.edu/documentation/nr/bookcpdf.html if required, im using the ran2 function.

The function should run fine, as its straight from the books data cd, but im not sure how i get it running.

It always returns the same (large) number. It says i initialize with a negative integer call to "idum" but im not sure how to do this.

Hope its all commented easy and i do the code tages right!

Code:

---

/* random number between 0 and 1 using ran2 method p282*/

#include <stdio.h>
#include <stdlib.h>

#define IM1 2147483563
#define IM2 2147483399
#define AM (1.0/IM1)
#define IMM1 (IM1-1)
#define IA1 40014
#define IA2 40692
#define IQ1 53668
#define IQ2 52774
#define IR1 12211
#define IR2 3791
#define NTAB 32
#define NDIV (1+IMM1/NTAB)
#define EPS 1.2e-7
#define RNMX (1.0-EPS)

/*Long period (>2e10^18) random number generator of L'Ecuyer with Bays-Durham shuffle
and added safeguards. Returns a random deviate between 0.0 and 1.0 (exlclusive of the endpoint values).
Call with idum  a negative integer to initialize;
thereafter do not alter idum between successive deviates in a sequence.
RNMX should approxiamte the largest floating value that is less than 1. */

float ran2(long *idum);

/*this bit i entered i doubt it is even close to being correct
i always get the same number, which isnt even in the range 0.0 to 1.0!*/

int main()
{

  float answer;
  long *pointer;
  long initial = -0.2;        /*it says im supposed to initialise with a call to idum, am i doing that here?*/
  pointer=&initial;


  answer=ran2(pointer);

  printf ("%d", answer);

  printf ("\nPress ENTER to continue.\n");

  fflush (stdin);
  (void) getchar();
  return (0);

}
/*end of bit i entered!*/

/*this is the ran2 function listed in the book, straight from the cd*/

float ran2(long *idum)
{
        int j;
        long k;
        static long idum2=123456789;
        static long iy=0;
        static long iv[NTAB];
        float temp;

        if (*idum <= 0)
  {
                if (-(*idum) < 1) *idum=1;
                else *idum = -(*idum);
                idum2=(*idum);
                for (j=NTAB+7;j>=0;j--)
      {
                        k=(*idum)/IQ1;
                        *idum=IA1*(*idum-k*IQ1)-k*IR1;
                        if (*idum < 0) *idum += IM1;
                        if (j < NTAB) iv[j] = *idum;
                }
                iy=iv[0];
        }
        k=(*idum)/IQ1;
        *idum=IA1*(*idum-k*IQ1)-k*IR1;
        if (*idum < 0) *idum += IM1;
        k=idum2/IQ2;
        idum2=IA2*(idum2-k*IQ2)-k*IR2;
        if (idum2 < 0) idum2 += IM2;
        j=iy/NDIV;
        iy=iv[j]-idum2;
        iv[j] = *idum;
        if (iy < 1) iy += IMM1;
        if ((temp=AM*iy) > RNMX) return RNMX;
        else return temp;
}

#undef IM1
#undef IM2
#undef AM
#undef IMM1
#undef IA1
#undef IA2
#undef IQ1
#undef IQ2
#undef IR1
#undef IR2
#undef NTAB
#undef NDIV
#undef EPS
#undef RNMX


/* (C) Copr. 1986-92 Numerical Recipes Software ]y. */

---

Any help much appreciated!!!

Code:

---

float answer;
printf ("%d", answer);

---

The %d is for ints; use %f.

Code:

---

long initial = -0.2; /* it says im supposed to initialise with a
                        call to idum, am i doing that here? */

---

The value -0.2 gets truncated to integral type.

Quote:

---

Originally Posted by Amoreldor

hey guys thanks for reading this.. i know that its a simple error on my part as i dont know what im doing..

It always returns the same (large) number. It says i initialize with a negative integer call to "idum" but im not sure how to do this.

Any help much appreciated!!!

---

What it says is:

Quote:

---

Call with idum a negative integer to initialize; thereafter, do not alter idum between successive deviates in a sequence.

---

So, you should declare a variable. Set the variable to a negative value and call ran2 with its address as the argument. Do not change the value again in your program (ran2 changes it for each call).

Here is an example that uses time(0) to get an initial value (so that the program gives different numbers each time it is run). I have printed out the value of *idum at each step to give a little visibility to the process, but you don't need it for anything.

Note that since the return value is a float, you need %f to print it:

Code:

---

#include <stdio.h>
#include <time.h>

int main()
{

  float answer;
  long *pointer;        /* don't really need this; just pass &idum each time */
  long idum = -time(0); /* get a value for initialization of ran2()          */
  int i;
  pointer=&idum;
  printf("*pointer = %ld\n", *pointer);

  for (i = 0; i < 10; i++) {
    answer=ran2(pointer);
    printf ("%2d: %f, *pointer = %ld\n", i, answer, *pointer);
  }

  return (0);
}

---

Note that this statement in your program sets *initial to 0:

Code:

---

  long initial = -0.2

---

<edit>: You need to #include <time.h> in order to use time().</edit>

Regards,

Dave

Brilliant, thanks a lot! Sorry was such a silly mistake!

The author of the book actually says "within the limits of its floating-point precision, ran2 provides perfect random numbers; a practical definition of perfect being ill pay $1000 to the first reader that convinces me otherwise....."

In light of that comment, I dont know if i should ask another question, but here goes..

No matter what negative value i choose to initialize with , ie -0.2, -0.4, -.565, i get the same answer out 0.285381. Is that because its the first in a "sequence of numbers generated from a first point"?

In which case i think i should be ok. If not i dont understand why i cant do this.

Sorry - i didnt refresh fast enough to see the second reply before responding!

Quote:

---

Originally Posted by Amoreldor

Brilliant, thanks a lot! Sorry was such a silly mistake!


No matter what negative value i choose to initialize with , ie -0.2, -0.4, -.565, i get the same answer out 0.285381. Is that because its the first in a "sequence of numbers generated from a first point"?

In which case i think i should be ok. If not i dont understand why i cant do this.

---

Must be a negative integer, since the argument is a long (long int, not a float).

You are always initializing with zero, since trying to assign a negative fraction to long variable results in truncation. (There is some funky code in ran2() that treats zero the same as -1. Man --- where did he get that from? I mean, really...)

Regards,

Dave

lol, thanks again- very helpful. Sorry i keep jumping the gun. Really am bad at this! Something tells me i need more coffee... (integer values)

Awsome, is all working perfectly. Am very grateful, many many thanks.

Amoreldor

Ok... again many thanks for all the help i have previously received, am now struggling with a new problem...

I want to use one of the numbers generated from the ran2 program to initialize a new program to give me a random number that lies within a Gaussian distribution between +1 and -1.

Ive had a bit of a bash at it, trying to follow what you guys had done to help me out with the other number set.

Its actually cool that the ran2 function gives me 9numbers at a time as im working in 8dimensions... (haha, dont ask- we dont know either- just geometry) and to get 9 gaussians from 9 inputs would be awsome... (cant believe i just wrote that)

I know i can do these things using very simple functions- but the "randomness" of what i get out is very influencial when testing the theory, so apologies for my over-complication of the number generation (so he says!)

Here is the function i am trying to use, as on disc, except using my ran2 answers instead of the ran1 as stated in the book:

Code:

---

/*Normal (Gaussian) deviate p289*/
#include <math.h>
float gasdev(long *idum)
{
        float ran1(long *idum);
        static int iset=0;
        static float gset;
        float fac,rsq,v1,v2;

        if  (iset == 0) {
                do {
                        v1=2.0*ran1(idum)-1.0;
                        v2=2.0*ran1(idum)-1.0;
                        rsq=v1*v1+v2*v2;
                } while (rsq >= 1.0 || rsq == 0.0);
                fac=sqrt(-2.0*log(rsq)/rsq);
                gset=v1*fac;
                iset=1;
                return v2*fac;
        } else {
                iset=0;
                return gset;
        }
}
/* (C) Copr. 1986-92 Numerical Recipes Software ]y.*/

---

here is my final ran2 code with the help of the board here:

Code:

---

/* random number between 0 and 1 using ran2 method p282*/

#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#include <time.h>

#define IM1 2147483563
#define IM2 2147483399
#define AM (1.0/IM1)
#define IMM1 (IM1-1)
#define IA1 40014
#define IA2 40692
#define IQ1 53668
#define IQ2 52774
#define IR1 12211
#define IR2 3791
#define NTAB 32
#define NDIV (1+IMM1/NTAB)
#define EPS 1.2e-7
#define RNMX (1.0-EPS)

/*Long period (>2e10^18) random number generator of L'Ecuyer with Bays-Durham shuffle
and added safeguards. Returns a random deviate between 0.0 and 1.0 (exlclusive of the endpoint values).
Call with idum  a negative integer to initialize;
thereafter do not alter idum between successive deviates in a sequence.
RNMX should approxiamte the largest floating value that is less than 1. */

float ran2(long *idum);

int main()
{

  float answer;
  long *pointer;        /* don't really need this; just pass &idum each time */
  long idum = -time(0); /* get a value for initialization of ran2()          */
  int i;
  pointer=&idum;
  printf("*pointer = %ld\n", *pointer);

  for (i = 0; i < 10; i++) {
    answer=ran2(pointer);
    printf ("%2d: %f, *pointer = %ld\n", i, answer, *pointer);
  }
  printf ("\nPress ENTER to continue.\n");

  fflush (stdin);
  (void) getchar();

  return (0);
}


/*this is the ran2 function listed in the book, straight from the cd*/

float ran2(long *idum)
{
        int j;
        long k;
        static long idum2=123456789;
        static long iy=0;
        static long iv[NTAB];
        float temp;

        if (*idum <= 0)
  {
                if (-(*idum) < 1) *idum=1;
                else *idum = -(*idum);
                idum2=(*idum);
                for (j=NTAB+7;j>=0;j--)
      {
                        k=(*idum)/IQ1;
                        *idum=IA1*(*idum-k*IQ1)-k*IR1;
                        if (*idum < 0) *idum += IM1;
                        if (j < NTAB) iv[j] = *idum;
                }
                iy=iv[0];
        }
        k=(*idum)/IQ1;
        *idum=IA1*(*idum-k*IQ1)-k*IR1;
        if (*idum < 0) *idum += IM1;
        k=idum2/IQ2;
        idum2=IA2*(idum2-k*IQ2)-k*IR2;
        if (idum2 < 0) idum2 += IM2;
        j=iy/NDIV;
        iy=iv[j]-idum2;
        iv[j] = *idum;
        if (iy < 1) iy += IMM1;
        if ((temp=AM*iy) > RNMX) return RNMX;
        else return temp;
}

#undef IM1
#undef IM2
#undef AM
#undef IMM1
#undef IA1
#undef IA2
#undef IQ1
#undef IQ2
#undef IR1
#undef IR2
#undef NTAB
#undef NDIV
#undef EPS
#undef RNMX

---

And here is my awful and messy attempt at patching the two together.... *cringes*

Code:

---

#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#include <time.h>

float gasdev(long *idum);

int main()
{
        float gaussiannumber;

  long *pointer;

  long idum = answer;

  int i;

  pointer=&idum;
  printf("*pointer = %1d\n", *pointer);

        for (i = 0; i < 10; i++)
  {
    gaussiannumber = ran2(pointer);
    printf ("%2d: %f, *pointer = %ld\n", i, gaussiannumber, *pointer);
  }
  printf ("\nPress ENTER to continue.\n");

  fflush (stdin);
  (void) getchar();

  return (0);

}

{

  float ran2(long *idum);
        static int iset=0;
        static float gaussiannumber;
        float fac,rsq,v1,v2;

        if  (iset == 0) {
                do {
                        v1=2.0*ran2(idum)-1.0;
                        v2=2.0*ran2(idum)-1.0;
                        rsq=v1*v1+v2*v2;
                } while (rsq >= 1.0 || rsq == 0.0);
                fac=sqrt(-2.0*log(rsq)/rsq);
                gaussiannumber=v1*fac;
                iset=1;
                return v2*fac;
        } else {
                iset=0;
                return gaussiannumber;
        }
}

---

Am sorry about posting so much code, just want to make everything clear.

Any help is much appreciated- im not very good at this logical stuff... or at C programming in general (have been doing it for 2wks! just got "thrown" at me by my professor)

*smallest violin*

Why did you rename gasdev() to ran2()?

gasdev() calls ran2(), so you have

Code:

---

float ran2(long *idum)
{
...
...
  v1=2.0*ran2(idum)-1.0;
...
...
}

---

So when you call this ran2() it calls itself recursively until the program stack overflows and does nothing.

Here's the plan:

Put in the code for version of ran2() that you have tested.
Put in the book's code for gasdev(), except use ran2() in place of ran1() if you want to.
Here's a suggestion for main() to test things:

Code:

---

int main()
{
  float gaussiannumber;

  long idum = -12345678; /* or use -time(0) to get different numbers each run */

  int i;

  for (i = 0; i < 10; i++) {
    gaussiannumber = gasdev(&idum);
    printf ("%2d: %9.6f\n", i, gaussiannumber);
  }

  return (0);
}

---

Regards,

Dave

Thanks very much, to be brutally honest i dont know why i did that, trying to copy the one you guys's helped me to get running. Thanks for the help, ill give it a good go!

Apologies for my random long posts!

Quote:

---

Originally Posted by Amoreldor

I want to use one of the numbers generated from the ran2 program to initialize a new program to give me a random number that lies within a Gaussian distribution between +1 and -1.

Its actually cool that the ran2 function gives me 9numbers at a time as im working in 8dimensions... (haha, dont ask- we dont know either- just geometry) and to get 9 gaussians from 9 inputs would be awsome... (cant believe i just wrote that)

---

This makes no sense to me: each call to ran2() gives a number from a uniform distribution between 0 and 1 (that is what it tries to do, anyhow, and the author says that numbers from this function satisfy all statistical tests for this).

gasdev() from the book calls ran1() (you call ran2() --- OK) several times to implement the Box-Muller method for generating numbers with normal distribution (Gaussian with mean=0 and standard deviation=1).

What does this have to do with these statements of yours?



Regards,

Dave

I want gasdev to call on the ran2 random number between 0and1 so that the gaussian number i can get out will be more "reliable" in its randomness.

In the problem im trying to model, we have to insert a one linear random number (0to1) and one gaussian number(-1to+1) so that we can start off a broyden method of root finding for our rather large and complicated free energy description of a spin glass. The idea is that we can then describe the turning points of the free energy with a hessain matrix for our solution.. phew..

Is that any clearer? As i said i dont understand the C-language and so not sure what im typing out- so im not suprised im confusing! But i really do appreciate the help!

p.s. i only really need one random linear and one gaussian to start things off, im just trying to be clever and close on the turning points from multiple starting values (hence the 9numbers comment)

Quote:

---

Originally Posted by Amoreldor

I want gasdev to call on the ran2 random number between 0and1 so that the gaussian number i can get out will be more "reliable" in its randomness.

In the problem im trying to model, we have to insert a one linear random number (0to1) and one gaussian number(-1to+1) so that we can start off a broyden method of root finding for our rather large and complicated free energy description of a spin glass. The idea is that we can then describe the turning points of the free energy with a hessain matrix for our solution.. phew..

Is that any clearer? As i said i dont understand the C-language and so not sure what im typing out- so im not suprised im confusing! But i really do appreciate the help!

p.s. i only really need one random linear and one gaussian to start things off, im just trying to be clever and close on the turning points from multiple starting values (hence the 9numbers comment)

---

My non-understanding of your statements has nothing to do with programming or programming language. I'm just trying to understand your terminology.

What do you mean by "gaussian number(-1 to 1)". gasdev() generates numbers from a normal (Gaussian) distribution. If I call gasdev(), say, 1000000 times, I get numbers between -4.6 and 4.6 (approximately). Sample mean is something like -.000859, sample variance is about 0.999. You can do other statisical tests.

Regards,

Dave

Sorry, i want to saddle the gaussian distribution about the f(x)=o vertical axis, limited between x=+1 and x=-1 (on a conventional x,y plot.)
Then i want to generate a random number that will lie upon that distribution. This is going to approximate my spin coupling which i will put into my free energy equation along with the ran2 generated number.

Apologies if im being unclear, im just starting the project so am still clarifying what im supposed to be researching too.. I have to run to circuit training now, so i wont be able to reply for an hour or so. Thanks for the help though!!