random_util.c

关于遗传算法的一些见地。特别是关于简单遗传程序设计的实现。

/**********************************************************************
  random_double_full()
  synopsis:	Return a random double within the allowed range.
  parameters:
  return:
  last updated:	16/06/01
 **********************************************************************/

double random_double_full(void)
  {
  return ( ((double)random_rand()/(double)RANDOM_RAND_MAX)*
                (DBL_MAX-DBL_MIN)+DBL_MIN );
  }


/**********************************************************************
  random_double()
  synopsis:	Return a random double within the specified range.
  parameters:
  return:
  last updated:	22/01/01
 **********************************************************************/

double random_double(const double max)
  {
  return ( max*(((double)random_rand())/(double)RANDOM_RAND_MAX) );
  }


/**********************************************************************
  random_double_range()
  synopsis:	Return a random double within the specified range.
  parameters:
  return:
  last updated:	07/06/00
 **********************************************************************/

double random_double_range(const double min, const double max)
  {
  return ( (max-min)*(((double)random_rand())/(double)RANDOM_RAND_MAX) + min );
  }


/**********************************************************************
  random_double_1()
  synopsis:	Return a random double within the range -1.0=>r>1.0
  parameters:
  return:
  last updated:	21/02/01
 **********************************************************************/

double random_double_1(void)
  {
  return ( 2.0*(((double)random_rand())/(double)RANDOM_RAND_MAX) - 1.0 );
  }


/**********************************************************************
  random_float_full()
  synopsis:	Return a random float within the allowed range.
  parameters:
  return:
  last updated:	4 Dec 2001
 **********************************************************************/

float random_float_full(void)
  {
  return ( ((float)random_rand()/(float)RANDOM_RAND_MAX)*
                (DBL_MAX-DBL_MIN)+DBL_MIN );
  }


/**********************************************************************
  random_float()
  synopsis:	Return a random float within the specified range.
  parameters:
  return:
  last updated:	4 Dec 2001
 **********************************************************************/

float random_float(const float max)
  {
  return ( max*(((float)random_rand())/(float)RANDOM_RAND_MAX) );
  }


/**********************************************************************
  random_float_range()
  synopsis:	Return a random float within the specified range.
  parameters:
  return:
  last updated:	4 Dec 2001
 **********************************************************************/

float random_float_range(const float min, const float max)
  {
  return ( (max-min)*(((float)random_rand())/(float)RANDOM_RAND_MAX) + min );
  }


/**********************************************************************
  random_float_1()
  synopsis:	Return a random float within the range -1.0=>r>1.0
  parameters:
  return:
  last updated:	4 Dec 2001
 **********************************************************************/

float random_float_1(void)
  {
  return ( 2.0*(((float)random_rand())/(float)RANDOM_RAND_MAX) - 1.0 );
  }


/**********************************************************************
  random_float_unit_uniform()
  synopsis:	Return a pseudo-random number with a uniform
		distribution in the range 0.0=>r>1.0
  parameters:
  return:
  last updated:	4 Dec 2001
 **********************************************************************/

float random_float_unit_uniform(void)
  {
  return ( (((float)random_rand())/(float)RANDOM_RAND_MAX) );
  }


/**********************************************************************
  random_unit_uniform()
  synopsis:	Return a pseudo-random number with a uniform
		distribution in the range 0.0=>r>1.0
  parameters:
  return:
  last updated:	11/01/01
 **********************************************************************/

double random_unit_uniform(void)
  {
  return ( (((double)random_rand())/(double)RANDOM_RAND_MAX) );
  }


/**********************************************************************
  random_float_gaussian()
  synopsis:	Return a pseudo-random number with a normal
		distribution with a given mean and standard devaiation.
  parameters:
  return:
  last updated:	4 Dec 2001
 **********************************************************************/

float random_float_gaussian(const float mean, const float stddev)
  {
  float	q,u,v,x,y;

/* 
 * Generate P = (u,v) uniform in rectangular acceptance region.
 */
  do
    {
    u = 1.0-random_float_unit_uniform();	/* in range 0.0>u>=1.0 */
    v = 1.7156 * (0.5 - random_float_unit_uniform());	/* in range 0.8578>v>=0.8578 */

/* Evaluate the quadratic form. */
    x = u - 0.449871;
    y = fabs(v) + 0.386595;
    q = x * x + y * (0.19600 * y - 0.25472 * x);

/*
 * Accept P if inside inner ellipse.
 * Reject P if outside outer ellipse, or outside acceptance region.
 */
    } while ((q >= 0.27597) && ((q > 0.27846) || (v * v > -4.0 * log(u) * u * u)));

/* Return ratio of P's coordinates as the normal deviate. */
  return (mean + 2.0 * stddev * v / u);	/* I'm not entirely sure why this *2.0 factor is needed! */
  }


/**********************************************************************
  random_float_unit_gaussian()
  synopsis:	Random number with normal distribution, average 0.0,
		deviation 1.0
  parameters:
  return:
  last updated:	07 Jan 2002
 **********************************************************************/

float random_float_unit_gaussian(void)
  {
  float			r, u, v, fac;
  static boolean	set = FALSE;
  static float		dset;

  if (set)
    {
    set = FALSE;
    return dset;
    }

  do
    {
    u = 2.0 * random_float_unit_uniform() - 1.0;
    v = 2.0 * random_float_unit_uniform() - 1.0;
    r = u*u + v*v;
    } while (r >= 1.0);

  fac = sqrt(-2.0 * log(r) / r);
  dset = v*fac;

  return u*fac;
  }

/* The original (thread-safe) version was this: */
#if 0
float random_float_unit_gaussian(void)
  {
  float	r, u, v;

  do
    {
    u = 2.0 * random_float_unit_uniform() - 1.0;
    v = 2.0 * random_float_unit_uniform() - 1.0;
    r = u*u + v*v;
    } while (r >= 1.0);

  return u*sqrt(-2.0 * log(r) / r);
  }
#endif


/**********************************************************************
  random_float_cauchy()
  synopsis:	Random number with a Cauchy/Lorentzian distribution.
  parameters:
  return:
  last updated:	4 Dec 2001
 **********************************************************************/

float random_float_cauchy(void)
  {
  return tan(random_float_range(-PI/2,PI/2));
  }


/**********************************************************************
  random_float_exponential()
  synopsis:	Random number with an exponential distribution, mean
		of 1.0.
  parameters:
  return:
  last updated:	4 Dec 2001
 **********************************************************************/

float random_float_exponential(void)
  {
  return -log(random_float_unit_uniform());
  }


/**********************************************************************
  random_gaussian()
  synopsis:	Return a pseudo-random number with a normal
		distribution with a given mean and standard devaiation.
  parameters:
  return:
  last updated:	11/01/01
 **********************************************************************/

/*
  Kinda based on: (But optimised quite a bit)

  ALGORITHM 712, COLLECTED ALGORITHMS FROM ACM.
  THIS WORK PUBLISHED IN TRANSACTIONS ON MATHEMATICAL SOFTWARE,
  VOL. 18, NO. 4, DECEMBER, 1992, PP. 434-435.
  The algorithm uses the ratio of uniforms method of A.J. Kinderman
  and J.F. Monahan augmented with quadratic bounding curves.
 */
double random_gaussian(const double mean, const double stddev)
  {
  double	q,u,v,x,y;

/* 
 * Generate P = (u,v) uniform in rectangular acceptance region.
 */
  do
    {
    u = 1.0-random_unit_uniform();	/* in range 0.0>u>=1.0 */
    v = 1.7156 * (0.5 - random_unit_uniform());	/* in range 0.8578>v>=0.8578 */

/* Evaluate the quadratic form. */
    x = u - 0.449871;
    y = fabs(v) + 0.386595;
    q = x * x + y * (0.19600 * y - 0.25472 * x);

/*
 * Accept P if inside inner ellipse.
 * Reject P if outside outer ellipse, or outside acceptance region.
 */
    } while ((q >= 0.27597) && ((q > 0.27846) || (v * v > -4.0 * log(u) * u * u)));

/* Return ratio of P's coordinates as the normal deviate. */
  return (mean + 2.0 * stddev * v / u);	/* I'm not entirely sure why this *2.0 factor is needed! */
  }


#if 0
/* Random number with normal distribution, average 0, deviation 1.
   From Numerical Recipes. */
double random_unit_gaussian(void)
  {
  static boolean	set = FALSE;
  static double		dset;
  double		fac, r, u, v;

  if (!set)
    {
    do
      {
      u = 2.0 * random_unit_uniform() - 1.0;
      v = 2.0 * random_unit_uniform() - 1.0;
      r = u*u + v*v;
      } while (r >= 1.0);

    fac = sqrt(-2.0 * log(r) / r);
    dset = u * fac;
    set = TRUE;

    return v * fac;
    }
  else
    {
    set = FALSE;
    return dset;
    }
  }
#endif


/**********************************************************************
  random_unit_gaussian()
  synopsis:	Random number with normal distribution, average 0.0,
		deviation 1.0
  parameters:
  return:
  last updated:	07 Jan 2002
 **********************************************************************/

double random_unit_gaussian(void)
  {
  double		r, u, v, fac;
  static boolean	set = FALSE;
  static double		dset;

  if (set)
    {
    set = FALSE;
    return dset;
    }

  do
    {
    u = 2.0 * random_unit_uniform() - 1.0;
    v = 2.0 * random_unit_uniform() - 1.0;
    r = u*u + v*v;
    } while (r >= 1.0);

  fac = sqrt(-2.0 * log(r) / r);
  dset = v*fac;

  return u*fac;
  }

/* The original (thread-safe) version was this: */
#if 0
double random_unit_gaussian(void)
  {
  double	r, u, v;

  do
    {
    u = 2.0 * random_unit_uniform() - 1.0;
    v = 2.0 * random_unit_uniform() - 1.0;
    r = u*u + v*v;
    } while (r >= 1.0);

  return u*sqrt(-2.0 * log(r) / r);
  }
#endif


/**********************************************************************
  random_cauchy()
  synopsis:	Random number with a Cauchy/Lorentzian distribution.
  parameters:	none
  return:	double	Random value.