|
7
|
1 // f-rand.cc -*- C++ -*- |
|
1
|
2 /* |
|
|
3 |
|
|
4 Copyright (C) 1993 John W. Eaton |
|
|
5 |
|
|
6 This file is part of Octave. |
|
|
7 |
|
|
8 Octave is free software; you can redistribute it and/or modify it |
|
|
9 under the terms of the GNU General Public License as published by the |
|
|
10 Free Software Foundation; either version 2, or (at your option) any |
|
|
11 later version. |
|
|
12 |
|
|
13 Octave is distributed in the hope that it will be useful, but WITHOUT |
|
|
14 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or |
|
|
15 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License |
|
|
16 for more details. |
|
|
17 |
|
|
18 You should have received a copy of the GNU General Public License |
|
|
19 along with Octave; see the file COPYING. If not, write to the Free |
|
|
20 Software Foundation, 675 Mass Ave, Cambridge, MA 02139, USA. |
|
|
21 |
|
|
22 */ |
|
|
23 |
|
|
24 #ifdef __GNUG__ |
|
|
25 #pragma implementation |
|
|
26 #endif |
|
|
27 |
|
|
28 #include "tree-const.h" |
|
|
29 #include "f77-uscore.h" |
|
|
30 #include "error.h" |
|
|
31 #include "utils.h" |
|
7
|
32 #include "f-rand.h" |
|
1
|
33 |
|
|
34 // Possible distributions of random numbers. |
|
|
35 enum rand_dist { uniform, normal }; |
|
|
36 |
|
|
37 // Current distribution of random numbers. |
|
|
38 static rand_dist current_distribution = uniform; |
|
|
39 |
|
|
40 extern "C" |
|
|
41 { |
|
|
42 int *F77_FCN (dgennor) (double*, double*, double*); |
|
|
43 int *F77_FCN (dgenunf) (double*, double*, double*); |
|
|
44 int *F77_FCN (setall) (int*, int*); |
|
|
45 int *F77_FCN (getsd) (int*, int*); |
|
|
46 } |
|
|
47 |
|
|
48 #ifdef WITH_DLD |
|
|
49 tree_constant * |
|
|
50 builtin_rand_2 (tree_constant *args, int nargin, int nargout) |
|
|
51 { |
|
|
52 return rand_internal (args, nargin, nargout); |
|
|
53 } |
|
|
54 #endif |
|
|
55 |
|
|
56 static double |
|
|
57 curr_rand_seed (void) |
|
|
58 { |
|
|
59 union d2i { double d; int i[2]; }; |
|
|
60 union d2i u; |
|
|
61 F77_FCN (getsd) (&(u.i[0]), &(u.i[1])); |
|
|
62 return u.d; |
|
|
63 } |
|
|
64 |
|
|
65 static int |
|
|
66 force_to_fit_range (int i, int lo, int hi) |
|
|
67 { |
|
|
68 assert (hi > lo && lo >= 0 && hi > lo); |
|
|
69 |
|
|
70 i = i > 0 ? i : -i; |
|
|
71 |
|
|
72 if (i < lo) |
|
|
73 i = lo; |
|
|
74 else if (i > hi) |
|
|
75 i = i % hi; |
|
|
76 |
|
|
77 return i; |
|
|
78 } |
|
|
79 |
|
|
80 static void |
|
|
81 set_rand_seed (double val) |
|
|
82 { |
|
|
83 union d2i { double d; int i[2]; }; |
|
|
84 union d2i u; |
|
|
85 u.d = val; |
|
|
86 int i0 = force_to_fit_range (u.i[0], 1, 2147483563); |
|
|
87 int i1 = force_to_fit_range (u.i[1], 1, 2147483399); |
|
|
88 F77_FCN (setall) (&i0, &i1); |
|
|
89 } |
|
|
90 |
|
|
91 static char * |
|
|
92 curr_rand_dist (void) |
|
|
93 { |
|
|
94 if (current_distribution == uniform) |
|
|
95 return "uniform"; |
|
|
96 else if (current_distribution == normal) |
|
|
97 return "normal"; |
|
|
98 else |
|
|
99 { |
|
|
100 panic_impossible (); |
|
|
101 return (char *) NULL; |
|
|
102 } |
|
|
103 } |
|
|
104 |
|
|
105 tree_constant * |
|
|
106 rand_internal (tree_constant *args, int nargin, int nargout) |
|
|
107 { |
|
|
108 // Assumes that we have been given the correct number of arguments. |
|
|
109 |
|
|
110 tree_constant *retval = NULL_TREE_CONST; |
|
|
111 |
|
|
112 static int initialized = 0; |
|
|
113 if (! initialized) |
|
|
114 { |
|
|
115 // Make the random number generator give us a different sequence every |
|
|
116 // time we start octave unless we specifically set the seed. |
|
|
117 #if 0 |
|
|
118 int s0 = 1234567890; |
|
|
119 int s1 = 123456789; |
|
|
120 #else |
|
|
121 time_t now; |
|
|
122 struct tm *tm; |
|
|
123 |
|
|
124 time (&now); |
|
|
125 tm = localtime (&now); |
|
|
126 |
|
|
127 int s0 = tm->tm_min * 60 + tm->tm_sec; |
|
|
128 int s1 = (tm->tm_mday - 1) * 24 * 3600 + tm->tm_hour * 3600 + s0; |
|
|
129 #endif |
|
|
130 s0 = force_to_fit_range (s0, 1, 2147483563); |
|
|
131 s1 = force_to_fit_range (s1, 1, 2147483399); |
|
|
132 |
|
|
133 F77_FCN (setall) (&s0, &s1); |
|
|
134 initialized = 1; |
|
|
135 } |
|
|
136 |
|
|
137 int n = 0; |
|
|
138 int m = 0; |
|
|
139 if (nargin == 1) |
|
|
140 { |
|
|
141 n = 1; |
|
|
142 m = 1; |
|
|
143 goto gen_matrix; |
|
|
144 } |
|
|
145 else if (nargin == 2) |
|
|
146 { |
|
|
147 switch (args[1].const_type ()) |
|
|
148 { |
|
|
149 case tree_constant_rep::string_constant: |
|
|
150 char *s_arg = args[1].string_value (); |
|
|
151 if (strcmp (s_arg, "dist") == 0) |
|
|
152 { |
|
|
153 retval = new tree_constant [2]; |
|
|
154 char *s = curr_rand_dist (); |
|
|
155 retval[0] = tree_constant (s); |
|
|
156 } |
|
|
157 else if (strcmp (s_arg, "seed") == 0) |
|
|
158 { |
|
|
159 retval = new tree_constant [2]; |
|
|
160 double d = curr_rand_seed (); |
|
|
161 retval[0] = tree_constant (d); |
|
|
162 } |
|
|
163 else if (strcmp (s_arg, "uniform") == 0) |
|
|
164 current_distribution = uniform; |
|
|
165 else if (strcmp (s_arg, "normal") == 0) |
|
|
166 current_distribution = normal; |
|
|
167 else |
|
|
168 { |
|
|
169 delete [] retval; |
|
|
170 retval = NULL_TREE_CONST; |
|
|
171 message ("rand", "unrecognized string argument"); |
|
|
172 } |
|
|
173 break; |
|
|
174 case tree_constant_rep::scalar_constant: |
|
|
175 case tree_constant_rep::complex_scalar_constant: |
|
|
176 n = NINT (args[1].double_value ()); |
|
|
177 m = n; |
|
|
178 goto gen_matrix; |
|
|
179 case tree_constant_rep::range_constant: |
|
|
180 { |
|
|
181 Range r = args[1].range_value (); |
|
|
182 n = 1; |
|
|
183 m = NINT (r.nelem ()); |
|
|
184 } |
|
|
185 goto gen_matrix; |
|
|
186 case tree_constant_rep::matrix_constant: |
|
|
187 case tree_constant_rep::complex_matrix_constant: |
|
|
188 n = NINT (args[1].rows ()); |
|
|
189 m = NINT (args[1].columns ()); |
|
|
190 goto gen_matrix; |
|
|
191 default: |
|
|
192 panic_impossible (); |
|
|
193 break; |
|
|
194 } |
|
|
195 } |
|
|
196 else if (nargin == 3) |
|
|
197 { |
|
|
198 if (args[1].is_string_type () |
|
|
199 && strcmp (args[1].string_value (), "seed") == 0) |
|
|
200 { |
|
|
201 double d = args[2].to_scalar (); |
|
|
202 set_rand_seed (d); |
|
|
203 } |
|
|
204 else |
|
|
205 { |
|
|
206 n = NINT (args[1].to_scalar ()); |
|
|
207 m = NINT (args[2].to_scalar ()); |
|
|
208 goto gen_matrix; |
|
|
209 } |
|
|
210 } |
|
|
211 |
|
|
212 return retval; |
|
|
213 |
|
|
214 gen_matrix: |
|
|
215 |
|
|
216 if (n == 0 || m == 0) |
|
|
217 { |
|
|
218 retval = new tree_constant [2]; |
|
|
219 Matrix m (0, 0); |
|
|
220 retval[0] = tree_constant (m); |
|
|
221 } |
|
|
222 else if (n > 0 && m > 0) |
|
|
223 { |
|
|
224 retval = new tree_constant [2]; |
|
|
225 Matrix rand_mat (n, m); |
|
|
226 for (int j = 0; j < m; j++) |
|
|
227 for (int i = 0; i < n; i++) |
|
|
228 { |
|
|
229 double d_zero = 0.0; |
|
|
230 double d_one = 1.0; |
|
|
231 double val; |
|
|
232 switch (current_distribution) |
|
|
233 { |
|
|
234 case uniform: |
|
|
235 F77_FCN (dgenunf) (&d_zero, &d_one, &val); |
|
|
236 rand_mat.elem (i, j) = val; |
|
|
237 break; |
|
|
238 case normal: |
|
|
239 F77_FCN (dgennor) (&d_zero, &d_one, &val); |
|
|
240 rand_mat.elem (i, j) = val; |
|
|
241 break; |
|
|
242 default: |
|
|
243 panic_impossible (); |
|
|
244 break; |
|
|
245 } |
|
|
246 } |
|
|
247 |
|
|
248 retval[0] = tree_constant (rand_mat); |
|
|
249 } |
|
|
250 else |
|
|
251 message ("rand", "invalid negative argument"); |
|
|
252 |
|
|
253 return retval; |
|
|
254 } |
|
|
255 |
|
|
256 /* |
|
|
257 ;;; Local Variables: *** |
|
|
258 ;;; mode: C++ *** |
|
|
259 ;;; page-delimiter: "^/\\*" *** |
|
|
260 ;;; End: *** |
|
|
261 */ |
