1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
|
SUBROUTINE ZLARNV( IDIST, ISEED, N, X )
*
* -- LAPACK auxiliary routine (version 3.1) --
* Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd..
* November 2006
*
* .. Scalar Arguments ..
INTEGER IDIST, N
* ..
* .. Array Arguments ..
INTEGER ISEED( 4 )
COMPLEX*16 X( * )
* ..
*
* Purpose
* =======
*
* ZLARNV returns a vector of n random complex numbers from a uniform or
* normal distribution.
*
* Arguments
* =========
*
* IDIST (input) INTEGER
* Specifies the distribution of the random numbers:
* = 1: real and imaginary parts each uniform (0,1)
* = 2: real and imaginary parts each uniform (-1,1)
* = 3: real and imaginary parts each normal (0,1)
* = 4: uniformly distributed on the disc abs(z) < 1
* = 5: uniformly distributed on the circle abs(z) = 1
*
* ISEED (input/output) INTEGER array, dimension (4)
* On entry, the seed of the random number generator; the array
* elements must be between 0 and 4095, and ISEED(4) must be
* odd.
* On exit, the seed is updated.
*
* N (input) INTEGER
* The number of random numbers to be generated.
*
* X (output) COMPLEX*16 array, dimension (N)
* The generated random numbers.
*
* Further Details
* ===============
*
* This routine calls the auxiliary routine DLARUV to generate random
* real numbers from a uniform (0,1) distribution, in batches of up to
* 128 using vectorisable code. The Box-Muller method is used to
* transform numbers from a uniform to a normal distribution.
*
* =====================================================================
*
* .. Parameters ..
DOUBLE PRECISION ZERO, ONE, TWO
PARAMETER ( ZERO = 0.0D+0, ONE = 1.0D+0, TWO = 2.0D+0 )
INTEGER LV
PARAMETER ( LV = 128 )
DOUBLE PRECISION TWOPI
PARAMETER ( TWOPI = 6.2831853071795864769252867663D+0 )
* ..
* .. Local Scalars ..
INTEGER I, IL, IV
* ..
* .. Local Arrays ..
DOUBLE PRECISION U( LV )
* ..
* .. Intrinsic Functions ..
INTRINSIC DCMPLX, EXP, LOG, MIN, SQRT
* ..
* .. External Subroutines ..
EXTERNAL DLARUV
* ..
* .. Executable Statements ..
*
DO 60 IV = 1, N, LV / 2
IL = MIN( LV / 2, N-IV+1 )
*
* Call DLARUV to generate 2*IL real numbers from a uniform (0,1)
* distribution (2*IL <= LV)
*
CALL DLARUV( ISEED, 2*IL, U )
*
IF( IDIST.EQ.1 ) THEN
*
* Copy generated numbers
*
DO 10 I = 1, IL
X( IV+I-1 ) = DCMPLX( U( 2*I-1 ), U( 2*I ) )
10 CONTINUE
ELSE IF( IDIST.EQ.2 ) THEN
*
* Convert generated numbers to uniform (-1,1) distribution
*
DO 20 I = 1, IL
X( IV+I-1 ) = DCMPLX( TWO*U( 2*I-1 )-ONE,
$ TWO*U( 2*I )-ONE )
20 CONTINUE
ELSE IF( IDIST.EQ.3 ) THEN
*
* Convert generated numbers to normal (0,1) distribution
*
DO 30 I = 1, IL
X( IV+I-1 ) = SQRT( -TWO*LOG( U( 2*I-1 ) ) )*
$ EXP( DCMPLX( ZERO, TWOPI*U( 2*I ) ) )
30 CONTINUE
ELSE IF( IDIST.EQ.4 ) THEN
*
* Convert generated numbers to complex numbers uniformly
* distributed on the unit disk
*
DO 40 I = 1, IL
X( IV+I-1 ) = SQRT( U( 2*I-1 ) )*
$ EXP( DCMPLX( ZERO, TWOPI*U( 2*I ) ) )
40 CONTINUE
ELSE IF( IDIST.EQ.5 ) THEN
*
* Convert generated numbers to complex numbers uniformly
* distributed on the unit circle
*
DO 50 I = 1, IL
X( IV+I-1 ) = EXP( DCMPLX( ZERO, TWOPI*U( 2*I ) ) )
50 CONTINUE
END IF
60 CONTINUE
RETURN
*
* End of ZLARNV
*
END
|
