LAPACK 3.12.0
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches
slaran.f
Go to the documentation of this file.
1*> \brief \b SLARAN
2*
3* =========== DOCUMENTATION ===========
4*
5* Online html documentation available at
6* http://www.netlib.org/lapack/explore-html/
7*
8* Definition:
9* ===========
10*
11* REAL FUNCTION SLARAN( ISEED )
12*
13* .. Array Arguments ..
14* INTEGER ISEED( 4 )
15* ..
16*
17*
18*> \par Purpose:
19* =============
20*>
21*> \verbatim
22*>
23*> SLARAN returns a random real number from a uniform (0,1)
24*> distribution.
25*> \endverbatim
26*
27* Arguments:
28* ==========
29*
30*> \param[in,out] ISEED
31*> \verbatim
32*> ISEED is INTEGER array, dimension (4)
33*> On entry, the seed of the random number generator; the array
34*> elements must be between 0 and 4095, and ISEED(4) must be
35*> odd.
36*> On exit, the seed is updated.
37*> \endverbatim
38*
39* Authors:
40* ========
41*
42*> \author Univ. of Tennessee
43*> \author Univ. of California Berkeley
44*> \author Univ. of Colorado Denver
45*> \author NAG Ltd.
46*
47*> \ingroup real_matgen
48*
49*> \par Further Details:
50* =====================
51*>
52*> \verbatim
53*>
54*> This routine uses a multiplicative congruential method with modulus
55*> 2**48 and multiplier 33952834046453 (see G.S.Fishman,
56*> 'Multiplicative congruential random number generators with modulus
57*> 2**b: an exhaustive analysis for b = 32 and a partial analysis for
58*> b = 48', Math. Comp. 189, pp 331-344, 1990).
59*>
60*> 48-bit integers are stored in 4 integer array elements with 12 bits
61*> per element. Hence the routine is portable across machines with
62*> integers of 32 bits or more.
63*> \endverbatim
64*>
65* =====================================================================
66 REAL function slaran( iseed )
67*
68* -- LAPACK auxiliary routine --
69* -- LAPACK is a software package provided by Univ. of Tennessee, --
70* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
71*
72* .. Array Arguments ..
73 INTEGER iseed( 4 )
74* ..
75*
76* =====================================================================
77*
78* .. Parameters ..
79 INTEGER m1, m2, m3, m4
80 parameter( m1 = 494, m2 = 322, m3 = 2508, m4 = 2549 )
81 REAL one
82 parameter( one = 1.0e+0 )
83 INTEGER ipw2
84 REAL r
85 parameter( ipw2 = 4096, r = one / ipw2 )
86* ..
87* .. Local Scalars ..
88 INTEGER it1, it2, it3, it4
89 REAL rndout
90* ..
91* .. Intrinsic Functions ..
92 INTRINSIC mod, real
93* ..
94* .. Executable Statements ..
95 10 CONTINUE
96*
97* multiply the seed by the multiplier modulo 2**48
98*
99 it4 = iseed( 4 )*m4
100 it3 = it4 / ipw2
101 it4 = it4 - ipw2*it3
102 it3 = it3 + iseed( 3 )*m4 + iseed( 4 )*m3
103 it2 = it3 / ipw2
104 it3 = it3 - ipw2*it2
105 it2 = it2 + iseed( 2 )*m4 + iseed( 3 )*m3 + iseed( 4 )*m2
106 it1 = it2 / ipw2
107 it2 = it2 - ipw2*it1
108 it1 = it1 + iseed( 1 )*m4 + iseed( 2 )*m3 + iseed( 3 )*m2 +
109 $ iseed( 4 )*m1
110 it1 = mod( it1, ipw2 )
111*
112* return updated seed
113*
114 iseed( 1 ) = it1
115 iseed( 2 ) = it2
116 iseed( 3 ) = it3
117 iseed( 4 ) = it4
118*
119* convert 48-bit integer to a real number in the interval (0,1)
120*
121 rndout = r*( real( it1 )+r*( real( it2 )+r*( real( it3 )+r*
122 $ ( real( it4 ) ) ) ) )
123*
124 IF (rndout.EQ.1.0) THEN
125* If a real number has n bits of precision, and the first
126* n bits of the 48-bit integer above happen to be all 1 (which
127* will occur about once every 2**n calls), then SLARAN will
128* be rounded to exactly 1.0. In IEEE single precision arithmetic,
129* this will happen relatively often since n = 24.
130* Since SLARAN is not supposed to return exactly 0.0 or 1.0
131* (and some callers of SLARAN, such as CLARND, depend on that),
132* the statistically correct thing to do in this situation is
133* simply to iterate again.
134* N.B. the case SLARAN = 0.0 should not be possible.
135*
136 GOTO 10
137 END IF
138*
139 slaran = rndout
140 RETURN
141*
142* End of SLARAN
143*
144 END
real function slaran(iseed)
SLARAN
Definition slaran.f:67