C: Random walk (Brownian motion) program does not return expected theoretical value. Why?

Question

This is just an experiment based on section 6-3 in "Feynman Lectures on Physics":

In its simplest version, we imagine a “game” in which a “player” starts at the point x=0 and at each “move” is required to take a step either forward (toward +x) or backward (toward −x). The choice is to be made randomly, determined, for example, by the toss of a coin.

Source: http://www.feynmanlectures.caltech.edu/I_06.html#Ch6-S3

My objective is to calculate the expected distance from the stating point. So, I suppose that each step is equal to one unit of distance. I wrote a simple C program to simulate 30 random steps, then calculate the final distance from the starting point. This is repeated for a million times, and the program averages the distance to get the expected distance.

Theoretically, the expected distance should be the square root of the number of the steps. That should be about sqrt(30) = 5.48.

However, the program is run few times and keeps returning a value near to 4.33 (to be more exact, 4.33461, 4.33453, and 4.34045). Why is it not even near to the theoretical value of about 5.48?

Here is my code:

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

int main ( int argc, char *argv[] )
{

  int number_of_steps = 30;
  int repetition = 1000000;
  int distance = 0;
  int total_distance = 0;
  double expected_distance;
  int i, j;

  srand(time(NULL));

  for ( i = 0; i < repetition; i++ ) {

    for ( j = 0; j < number_of_steps; j++) {
      distance += rand() & 1 ? -1 : 1;
    }

    total_distance += abs(distance);
    distance = 0;

  }

  expected_distance = (float) total_distance / i;

  printf ( "%g\n", expected_distance );
  return EXIT_SUCCESS;
}       /* ----------  end of function main  ---------- */

Answer 1

From the lecture you linked to, your theoretical expectation is based on the root mean square, which is different from the arithmetic mean, which is what you have coded. By changing the algorithm from one to the other, the code now gives you the expected results.

for ( i = 0; i < repetition; i++ ) {

    for ( j = 0; j < number_of_steps; j++) {
      distance += rand() & 1 ? -1 : 1;
    }

    total_distance += distance * distance;
    distance = 0;

  }

  expected_distance = sqrt((float) total_distance / repetition);

  printf ( "%g\n", expected_distance );
  return EXIT_SUCCESS;
}

Answer 2

The answer to this post suggests that using the low-order bit(s) of rand() is unlikely to be a great choice.

I'd try a different way of generating your +1 or -1.