Random numbers: Difference between revisions

From Rosetta Code
Content added Content deleted
(added Fortran)
m (Ada code cleanup)
Line 6: Line 6:


=={{header|Ada}}==
=={{header|Ada}}==
<Ada>
with Ada.Numerics.Float_Random; use Ada.Numerics.Float_Random;
with Ada.Numerics.Generic_Elementary_Functions;
with Ada.Numerics; use Ada.Numerics;
with Ada.Numerics.Float_Random; use Ada.Numerics.Float_Random;
with Ada.Numerics.Elementary_Functions; use Ada.Numerics.Elementary_Functions;
procedure Normal_Random is
Seed : Generator;
function Normal_Distribution(Seed : Generator) return Float is
package Elementary_Flt is new
Ada.Numerics.Generic_Elementary_Functions(Float);
use Elementary_Flt;
use Ada.Numerics;
R1 : Float;
R2 : Float;
Mu : constant Float := 1.0;
Sigma : constant Float := 0.5;
begin
R1 := Random(Seed);
R2 := Random(Seed);
return Mu + (Sigma * Sqrt(-2.0 * Log(X => R1, Base => 10.0)) *
Cos(2.0 * Pi * R2));
end Normal_Distribution;
type Normal_Array is array(1..1000) of Float;
Distribution : Normal_Array;
begin
Reset(Seed);
for I in Distribution'range loop
Distribution(I) := Normal_Distribution(Seed);
end loop;
end Normal_Random;


procedure Normal_Random is
function Normal_Distribution
( Seed : Generator;
Mu : Float := 1.0;
Sigma : Float := 0.5
) return Float is
begin
return
Mu + (Sigma * Sqrt (-2.0 * Log (Random (Seed), 10.0)) * Cos (2.0 * Pi * Random (Seed)));
end Normal_Distribution;
Seed : Generator;
Distribution : array (1..1_000) of Float;
begin
Reset (Seed);
for I in Distribution'Range loop
Distribution (I) := Normal_Distribution (Seed);
end loop;
end Normal_Random;
</Ada>
=={{header|BASIC}}==
=={{header|BASIC}}==
{{works with|QuickBasic|4.5}}
{{works with|QuickBasic|4.5}}

Revision as of 13:10, 7 June 2008

Task
Random numbers
You are encouraged to solve this task according to the task description, using any language you may know.

The goal of this task is to generate a collection filled with 1000 normally distributed random numbers with a mean of 1.0 and a standard deviation of 0.5

Many libraries only generate uniformly distributed random numbers. If so, use this formula to convert them to a normal distribution.

Ada

<Ada> with Ada.Numerics; use Ada.Numerics; with Ada.Numerics.Float_Random; use Ada.Numerics.Float_Random; with Ada.Numerics.Elementary_Functions; use Ada.Numerics.Elementary_Functions;

procedure Normal_Random is

  function Normal_Distribution
           (  Seed  : Generator;
              Mu    : Float := 1.0;
              Sigma : Float := 0.5
           )  return Float is 
  begin
     return
        Mu + (Sigma * Sqrt (-2.0 * Log (Random (Seed), 10.0)) * Cos (2.0 * Pi * Random (Seed)));
  end Normal_Distribution;
     
  Seed         : Generator;
  Distribution : array (1..1_000) of Float; 

begin

  Reset (Seed);
  for I in Distribution'Range loop
     Distribution (I) := Normal_Distribution (Seed);
  end loop;

end Normal_Random; </Ada>

BASIC

Works with: QuickBasic version 4.5
RANDOMIZE TIMER 'seeds random number generator with the system time
pi = 3.141592653589793#
DIM a(1 TO 1000) AS DOUBLE
CLS
FOR i = 1 TO 1000
   a(i) = 1 + SQR(-2 * LOG(RND)) * COS(2 * pi * RND)
NEXT i

C

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

double drand()   /* uniform distribution, (0..1] */
{
  return (rand()+1.0)/(RAND_MAX+1.0);
}
double random_normal()  /* normal distribution, centered on 0, std dev 1 */
{
  return sqrt(-2*log(drand())) * cos(2*M_PI*drand());
}
int main()
{
  int i;
  double rands[1000];
  for (i=0; i<1000; i++)
    rands[i] = 1.0 + 0.5*random_normal();
  return 0;
}

C++

#include <cstdlib>   // for rand
#include <cmath>     // for atan, sqrt, log, cos
#include <algorithm> // for generate_n

double const pi = 4*std::atan(1.0);

// simple functor for normal distribution
class normal_distribution
{
public:
  normal_distribution(double m, double s): mu(m), sigma(s) {}
  double operator() // returns a single normally distributed number
  {
    double r1 = (std::rand() + 1.0)/(RAND_MAX + 1.0); // gives equal distribution in (0, 1]
    double r2 = (std::rand() + 1.0)/(RAND_MAX + 1.0);
    return mu + sigma * std::sqrt(-2*std::log(r1))*std::cos(2*pi*r2);
  }
private:
  double mu, sigma;
};

int main()
{
  double array[1000];
  std::generate_n(array, 1000, normal_distribution(1.0, 0.5));
}

E

accum [] for _ in 1..1000 { _.with(entropy.nextGaussian()) }

Forth

Works with: gforth version 0.6.2
require random.fs
here to seed

-1. 1 rshift 2constant MAX-D	\ or s" MAX-D" ENVIRONMENT? drop

: frnd ( -- f )			\ uniform distribution 0..1
  rnd rnd dabs d>f MAX-D d>f f/ ;

: frnd-normal ( -- f )		\ centered on 0, std dev 1
  frnd pi f* 2e f* fcos
  frnd fln -2e f* fsqrt f* ;

: ,normals ( n -- )		\ store many, centered on 1, std dev 0.5
  0 do frnd-normal 0.5e f* 1e f+ f, loop ;

create rnd-array 1000 ,normals

Fortran

Works with: Fortran version 90 and later
PROGRAM Random

  INTEGER, PARAMETER :: n = 1000
  INTEGER :: i
  REAL :: array(n), pi, temp, mean = 1.0, sd = 0.5

  pi = 4.0*ATAN(1.0)
  CALL RANDOM_NUMBER(array) ! Uniform distribution
 
! Now convert to normal distribution
  DO i = 1, n-1, 2
    temp = sd * SQRT(-2.0*LOG(array(i))) * COS(2*pi*array(i+1)) + mean
    array(i+1) = sd * SQRT(-2.0*LOG(array(i))) * SIN(2*pi*array(i+1)) + mean
    array(i) = temp
  END DO

! Check mean and standard deviation
  mean = SUM(array)/n
  temp = 0.0
  DO i = 1, n
    temp = temp + (array(i) - mean)**2
  END DO
  sd = SQRT(temp/n)
 
  WRITE(*, "(A,F8.6)") "Mean = ", mean
  WRITE(*, "(A,F8.6)") "Standard Deviation = ", sd

END PROGRAM Random

Output

Mean = 0.995112
Standard Deviation = 0.503373

Haskell

import System.Random

pairs :: [a] -> [(a,a)]
pairs (x:y:zs) = (x,y):pairs zs
pairs _        = []

gauss mu sigma (r1,r2) = 
  mu + sigma * sqrt (-2 * log r1) * cos (2 * pi * r2)

gaussians :: (RandomGen g, Random a, Floating a) => Int -> g -> [a]
gaussians n g = take n $ map (gauss 1.0 0.5) $ pairs $ randoms g
result :: IO [Double]
result = getStdGen >>= \g -> return $ gaussians 1000 g

Groovy

rnd = new Random()
result = (1..1000).inject([]) { r, i -> r << rnd.nextGaussian() }

IDL

result = 1.0 + 0.5*randomn(seed,1000)

J

urand=: ?@$ 0: 
zrand=: (2 o. 2p1 * urand) * [: %: _2 * [: ^. urand

1 + 0.5 * zrand 100

Java

double[] list = new double[1000];
Random rng = new Random();
for(int i = 0;i<list.length;i++) {
  list[i] = 1.0 + 0.5 * rng.nextGaussian()
}

JavaScript

function randomNormal() {
  return Math.cos(2 * Math.PI * Math.random()) * Math.sqrt(-2 * Math.log(Math.random()));
}

var a = new Array(1000);
for (var i=0; i<a.length; i++)
  a[i] = randomNormal() / 2 + 1;

Works with: UCB Logo

The earliest Logos only have a RANDOM function for picking a random non-negative integer. Many modern Logos have floating point random generators built-in.

to random.float   ; 0..1
  localmake "max.int lshift -1 -1
  output quotient random :max.int :max.int
end

to random.gaussian
  output product cos random 360  sqrt -2 / ln random.float
end

make "randoms cascade 1000 [fput random.gaussian / 2 + 1 ?] []

MAXScript

arr = #()
for i in 1 to 1000 do 
(
    a = random 0.0 1.0
    b = random 0.0 1.0
    c = 1.0 + 0.5 * sqrt (-2*log a) * cos (360*b) -- Maxscript cos takes degrees
    append arr c
)

OCaml

let pi = 4. *. atan 1.;;
let random_gaussian () =
  1. +. sqrt (-2. *. log (Random.float 1.)) *. cos (2. *. pi *. Random.float 1.);;
let a = Array.init 1000 (fun _ -> random_gaussian ());;

Perl

use Math::Cephes qw($PI);

map {
    1.0 + 0.5 * sqrt (-2 * log rand) * cos (2 * $PI * rand)
} 1..1000

PHP

$pi 	= pi();          // Set PI
$a	= range(1,1000); // Create array
// Cycle array values
foreach(range(1,1000) as $i){
      $a[$i] = 1 + sqrt(-2 * log(mt_rand())) * cos(2 * $pi * mt_rand());
}

Pop11

;;; Choose radians as arguments to trigonometic functions
true -> popradians;
;;; procedure generating standard normal distribution
define random_normal() -> result;
lvars r1 = random0(1.0), r2 = random0(1.0);
     cos(2*pi*r1)*sqrt(-2*log(r2)) -> result
enddefine;
lvars array, i;
;;; Put numbers on the stack
for i from 1 to 1000 do 1.0+0.5*random_normal() endfor;
;;; collect them into array
consvector(1000) -> array;

Python

Works with: Python version 2.5
import random
randList = [random.gauss(1, .5) for i in range(1000)]
# or [ random.normalvariate(1, 0.5) for i in range(1000)]

Note that the random module in the Python standard library supports a number of statistical distribution methods.

Tcl

proc nrand {} {return [expr sqrt(-2*log(rand()))*cos(4*acos(0)*rand())]}
for {set i 0} {$i < 1000} {incr i} {lappend result [expr 1+.5*nrand()]}

TI-83 BASIC

Calculator symbol translations:

"STO" arrow: →

Square root sign: √

ClrList L1
Radian
For(A,1,1000)
√(-2*ln(rand))*cos(2*π*A)→L1(A)
End