1 /* ./src_f77/clarnv.f -- translated by f2c (version 20030320).
2 You must link the resulting object file with the libraries:
3 -lf2c -lm (in that order)
4 */
5
6 #include <punc/vf2c.h>
7
clarnv_(integer * idist,integer * iseed,integer * n,complex * x)8 /* Subroutine */ int clarnv_(integer *idist, integer *iseed, integer *n,
9 complex *x)
10 {
11 /* System generated locals */
12 integer i__1, i__2, i__3, i__4, i__5;
13 real r__1, r__2;
14 complex q__1, q__2, q__3;
15
16 /* Builtin functions */
17 double log(doublereal), sqrt(doublereal);
18 void c_exp(complex *, complex *);
19
20 /* Local variables */
21 static integer i__;
22 static real u[128];
23 static integer il, iv;
24 extern /* Subroutine */ int slaruv_(integer *, integer *, real *);
25
26
27 /* -- LAPACK auxiliary routine (version 3.0) -- */
28 /* Univ. of Tennessee, Univ. of California Berkeley, NAG Ltd., */
29 /* Courant Institute, Argonne National Lab, and Rice University */
30 /* September 30, 1994 */
31
32 /* .. Scalar Arguments .. */
33 /* .. */
34 /* .. Array Arguments .. */
35 /* .. */
36
37 /* Purpose */
38 /* ======= */
39
40 /* CLARNV returns a vector of n random complex numbers from a uniform or */
41 /* normal distribution. */
42
43 /* Arguments */
44 /* ========= */
45
46 /* IDIST (input) INTEGER */
47 /* Specifies the distribution of the random numbers: */
48 /* = 1: real and imaginary parts each uniform (0,1) */
49 /* = 2: real and imaginary parts each uniform (-1,1) */
50 /* = 3: real and imaginary parts each normal (0,1) */
51 /* = 4: uniformly distributed on the disc abs(z) < 1 */
52 /* = 5: uniformly distributed on the circle abs(z) = 1 */
53
54 /* ISEED (input/output) INTEGER array, dimension (4) */
55 /* On entry, the seed of the random number generator; the array */
56 /* elements must be between 0 and 4095, and ISEED(4) must be */
57 /* odd. */
58 /* On exit, the seed is updated. */
59
60 /* N (input) INTEGER */
61 /* The number of random numbers to be generated. */
62
63 /* X (output) COMPLEX array, dimension (N) */
64 /* The generated random numbers. */
65
66 /* Further Details */
67 /* =============== */
68
69 /* This routine calls the auxiliary routine SLARUV to generate random */
70 /* real numbers from a uniform (0,1) distribution, in batches of up to */
71 /* 128 using vectorisable code. The Box-Muller method is used to */
72 /* transform numbers from a uniform to a normal distribution. */
73
74 /* ===================================================================== */
75
76 /* .. Parameters .. */
77 /* .. */
78 /* .. Local Scalars .. */
79 /* .. */
80 /* .. Local Arrays .. */
81 /* .. */
82 /* .. Intrinsic Functions .. */
83 /* .. */
84 /* .. External Subroutines .. */
85 /* .. */
86 /* .. Executable Statements .. */
87
88 /* Parameter adjustments */
89 --x;
90 --iseed;
91
92 /* Function Body */
93 i__1 = *n;
94 for (iv = 1; iv <= i__1; iv += 64) {
95 /* Computing MIN */
96 i__2 = 64, i__3 = *n - iv + 1;
97 il = min(i__2,i__3);
98
99 /* Call SLARUV to generate 2*IL real numbers from a uniform (0,1) */
100 /* distribution (2*IL <= LV) */
101
102 i__2 = il << 1;
103 slaruv_(&iseed[1], &i__2, u);
104
105 if (*idist == 1) {
106
107 /* Copy generated numbers */
108
109 i__2 = il;
110 for (i__ = 1; i__ <= i__2; ++i__) {
111 i__3 = iv + i__ - 1;
112 i__4 = (i__ << 1) - 2;
113 i__5 = (i__ << 1) - 1;
114 q__1.r = u[i__4], q__1.i = u[i__5];
115 x[i__3].r = q__1.r, x[i__3].i = q__1.i;
116 /* L10: */
117 }
118 } else if (*idist == 2) {
119
120 /* Convert generated numbers to uniform (-1,1) distribution */
121
122 i__2 = il;
123 for (i__ = 1; i__ <= i__2; ++i__) {
124 i__3 = iv + i__ - 1;
125 r__1 = u[(i__ << 1) - 2] * 2.f - 1.f;
126 r__2 = u[(i__ << 1) - 1] * 2.f - 1.f;
127 q__1.r = r__1, q__1.i = r__2;
128 x[i__3].r = q__1.r, x[i__3].i = q__1.i;
129 /* L20: */
130 }
131 } else if (*idist == 3) {
132
133 /* Convert generated numbers to normal (0,1) distribution */
134
135 i__2 = il;
136 for (i__ = 1; i__ <= i__2; ++i__) {
137 i__3 = iv + i__ - 1;
138 r__1 = sqrt(log(u[(i__ << 1) - 2]) * -2.f);
139 r__2 = u[(i__ << 1) - 1] * 6.2831853071795864769252867663f;
140 q__3.r = 0.f, q__3.i = r__2;
141 c_exp(&q__2, &q__3);
142 q__1.r = r__1 * q__2.r, q__1.i = r__1 * q__2.i;
143 x[i__3].r = q__1.r, x[i__3].i = q__1.i;
144 /* L30: */
145 }
146 } else if (*idist == 4) {
147
148 /* Convert generated numbers to complex numbers uniformly */
149 /* distributed on the unit disk */
150
151 i__2 = il;
152 for (i__ = 1; i__ <= i__2; ++i__) {
153 i__3 = iv + i__ - 1;
154 r__1 = sqrt(u[(i__ << 1) - 2]);
155 r__2 = u[(i__ << 1) - 1] * 6.2831853071795864769252867663f;
156 q__3.r = 0.f, q__3.i = r__2;
157 c_exp(&q__2, &q__3);
158 q__1.r = r__1 * q__2.r, q__1.i = r__1 * q__2.i;
159 x[i__3].r = q__1.r, x[i__3].i = q__1.i;
160 /* L40: */
161 }
162 } else if (*idist == 5) {
163
164 /* Convert generated numbers to complex numbers uniformly */
165 /* distributed on the unit circle */
166
167 i__2 = il;
168 for (i__ = 1; i__ <= i__2; ++i__) {
169 i__3 = iv + i__ - 1;
170 r__1 = u[(i__ << 1) - 1] * 6.2831853071795864769252867663f;
171 q__2.r = 0.f, q__2.i = r__1;
172 c_exp(&q__1, &q__2);
173 x[i__3].r = q__1.r, x[i__3].i = q__1.i;
174 /* L50: */
175 }
176 }
177 /* L60: */
178 }
179 return 0;
180
181 /* End of CLARNV */
182
183 } /* clarnv_ */
184
185