|
1
|
1 // Some extra friends of the tree constant class. -*- C++ -*- |
|
|
2 // See also the other tc-*.cc files. |
|
|
3 /* |
|
|
4 |
|
296
|
5 Copyright (C) 1992, 1993, 1994 John W. Eaton |
|
1
|
6 |
|
|
7 This file is part of Octave. |
|
|
8 |
|
|
9 Octave is free software; you can redistribute it and/or modify it |
|
|
10 under the terms of the GNU General Public License as published by the |
|
|
11 Free Software Foundation; either version 2, or (at your option) any |
|
|
12 later version. |
|
|
13 |
|
|
14 Octave is distributed in the hope that it will be useful, but WITHOUT |
|
|
15 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or |
|
|
16 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License |
|
|
17 for more details. |
|
|
18 |
|
|
19 You should have received a copy of the GNU General Public License |
|
|
20 along with Octave; see the file COPYING. If not, write to the Free |
|
|
21 Software Foundation, 675 Mass Ave, Cambridge, MA 02139, USA. |
|
|
22 |
|
|
23 */ |
|
|
24 |
|
240
|
25 #ifdef HAVE_CONFIG_H |
|
|
26 #include "config.h" |
|
1
|
27 #endif |
|
|
28 |
|
|
29 #include <strstream.h> |
|
|
30 #include <iostream.h> |
|
|
31 #include <fstream.h> |
|
|
32 #include <string.h> |
|
|
33 #include <ctype.h> |
|
|
34 |
|
456
|
35 #include "EIG.h" |
|
|
36 |
|
1
|
37 #include "unwind-prot.h" |
|
|
38 #include "tree-const.h" |
|
|
39 #include "user-prefs.h" |
|
|
40 #include "variables.h" |
|
|
41 #include "octave.h" |
|
|
42 #include "gripes.h" |
|
|
43 #include "error.h" |
|
|
44 #include "input.h" |
|
|
45 #include "octave-hist.h" |
|
|
46 #include "pager.h" |
|
|
47 #include "utils.h" |
|
|
48 #include "parse.h" |
|
|
49 #include "lex.h" |
|
|
50 |
|
|
51 Matrix |
|
|
52 max (const Matrix& a, const Matrix& b) |
|
|
53 { |
|
|
54 int nr = a.rows (); |
|
|
55 int nc = a.columns (); |
|
|
56 if (nr != b.rows () || nc != b.columns ()) |
|
226
|
57 { |
|
|
58 error ("two-arg max expecting args of same size"); |
|
|
59 return Matrix (); |
|
|
60 } |
|
1
|
61 |
|
|
62 Matrix result (nr, nc); |
|
|
63 |
|
|
64 for (int j = 0; j < nc; j++) |
|
|
65 for (int i = 0; i < nr; i++) |
|
|
66 { |
|
|
67 double a_elem = a.elem (i, j); |
|
|
68 double b_elem = b.elem (i, j); |
|
|
69 result.elem (i, j) = MAX (a_elem, b_elem); |
|
|
70 } |
|
|
71 |
|
|
72 return result; |
|
|
73 } |
|
|
74 |
|
199
|
75 ComplexMatrix |
|
|
76 max (const ComplexMatrix& a, const ComplexMatrix& b) |
|
|
77 { |
|
|
78 int nr = a.rows (); |
|
|
79 int nc = a.columns (); |
|
|
80 if (nr != b.rows () || nc != b.columns ()) |
|
226
|
81 { |
|
|
82 error ("two-arg max expecting args of same size"); |
|
|
83 return ComplexMatrix (); |
|
|
84 } |
|
199
|
85 |
|
|
86 ComplexMatrix result (nr, nc); |
|
|
87 |
|
|
88 for (int j = 0; j < nc; j++) |
|
|
89 for (int i = 0; i < nr; i++) |
|
|
90 { |
|
|
91 double abs_a_elem = abs (a.elem (i, j)); |
|
|
92 double abs_b_elem = abs (b.elem (i, j)); |
|
|
93 if (abs_a_elem > abs_b_elem) |
|
|
94 result.elem (i, j) = a.elem (i, j); |
|
|
95 else |
|
|
96 result.elem (i, j) = b.elem (i, j); |
|
|
97 } |
|
|
98 |
|
|
99 return result; |
|
|
100 } |
|
|
101 |
|
1
|
102 Matrix |
|
|
103 min (const Matrix& a, const Matrix& b) |
|
|
104 { |
|
|
105 int nr = a.rows (); |
|
|
106 int nc = a.columns (); |
|
|
107 if (nr != b.rows () || nc != b.columns ()) |
|
226
|
108 { |
|
|
109 error ("two-arg min expecting args of same size"); |
|
|
110 return Matrix (); |
|
|
111 } |
|
1
|
112 |
|
|
113 Matrix result (nr, nc); |
|
|
114 |
|
|
115 for (int j = 0; j < nc; j++) |
|
|
116 for (int i = 0; i < nr; i++) |
|
|
117 { |
|
|
118 double a_elem = a.elem (i, j); |
|
|
119 double b_elem = b.elem (i, j); |
|
|
120 result.elem (i, j) = MIN (a_elem, b_elem); |
|
|
121 } |
|
|
122 |
|
|
123 return result; |
|
|
124 } |
|
|
125 |
|
199
|
126 ComplexMatrix |
|
|
127 min (const ComplexMatrix& a, const ComplexMatrix& b) |
|
|
128 { |
|
|
129 int nr = a.rows (); |
|
|
130 int nc = a.columns (); |
|
|
131 if (nr != b.rows () || nc != b.columns ()) |
|
226
|
132 { |
|
|
133 error ("two-arg min expecting args of same size"); |
|
|
134 return ComplexMatrix (); |
|
|
135 } |
|
199
|
136 |
|
|
137 ComplexMatrix result (nr, nc); |
|
|
138 |
|
|
139 for (int j = 0; j < nc; j++) |
|
|
140 for (int i = 0; i < nr; i++) |
|
|
141 { |
|
|
142 double abs_a_elem = abs (a.elem (i, j)); |
|
|
143 double abs_b_elem = abs (b.elem (i, j)); |
|
|
144 if (abs_a_elem < abs_b_elem) |
|
|
145 result.elem (i, j) = a.elem (i, j); |
|
|
146 else |
|
|
147 result.elem (i, j) = b.elem (i, j); |
|
|
148 } |
|
|
149 |
|
|
150 return result; |
|
|
151 } |
|
|
152 |
|
1
|
153 static void |
|
164
|
154 get_dimensions (const tree_constant& a, const char *warn_for, |
|
|
155 int& nr, int& nc) |
|
1
|
156 { |
|
|
157 tree_constant tmpa = a.make_numeric (); |
|
|
158 |
|
|
159 if (tmpa.is_scalar_type ()) |
|
|
160 { |
|
|
161 double tmp = tmpa.double_value (); |
|
|
162 nr = nc = NINT (tmp); |
|
|
163 } |
|
|
164 else |
|
|
165 { |
|
481
|
166 nr = tmpa.rows (); |
|
|
167 nc = tmpa.columns (); |
|
|
168 |
|
|
169 if ((nr == 1 && nc == 2) || (nr == 2 && nc == 1)) |
|
|
170 { |
|
|
171 ColumnVector v = tmpa.to_vector (); |
|
|
172 |
|
|
173 nr = NINT (v.elem (0)); |
|
|
174 nc = NINT (v.elem (1)); |
|
|
175 } |
|
|
176 else |
|
|
177 warning ("%s (A): use %s (size (A)) instead", warn_for, warn_for); |
|
1
|
178 } |
|
|
179 |
|
143
|
180 check_dimensions (nr, nc, warn_for); // May set error_state. |
|
1
|
181 } |
|
|
182 |
|
|
183 static void |
|
164
|
184 get_dimensions (const tree_constant& a, const tree_constant& b, |
|
|
185 const char *warn_for, int& nr, int& nc) |
|
1
|
186 { |
|
|
187 tree_constant tmpa = a.make_numeric (); |
|
|
188 tree_constant tmpb = b.make_numeric (); |
|
|
189 |
|
|
190 if (tmpa.is_scalar_type () && tmpb.is_scalar_type ()) |
|
|
191 { |
|
|
192 nr = NINT (tmpa.double_value ()); |
|
|
193 nc = NINT (tmpb.double_value ()); |
|
|
194 |
|
143
|
195 check_dimensions (nr, nc, warn_for); // May set error_state. |
|
1
|
196 } |
|
|
197 else |
|
143
|
198 error ("%s: expecting two scalar arguments", warn_for); |
|
1
|
199 } |
|
|
200 |
|
|
201 tree_constant |
|
164
|
202 fill_matrix (const tree_constant& a, double val, const char *warn_for) |
|
1
|
203 { |
|
|
204 int nr, nc; |
|
143
|
205 get_dimensions (a, warn_for, nr, nc); |
|
|
206 |
|
|
207 if (error_state) |
|
|
208 return tree_constant (); |
|
1
|
209 |
|
|
210 Matrix m (nr, nc, val); |
|
|
211 |
|
|
212 return tree_constant (m); |
|
|
213 } |
|
|
214 |
|
|
215 tree_constant |
|
164
|
216 fill_matrix (const tree_constant& a, const tree_constant& b, |
|
|
217 double val, const char *warn_for) |
|
1
|
218 { |
|
|
219 int nr, nc; |
|
143
|
220 get_dimensions (a, b, warn_for, nr, nc); // May set error_state. |
|
|
221 |
|
|
222 if (error_state) |
|
|
223 return tree_constant (); |
|
1
|
224 |
|
|
225 Matrix m (nr, nc, val); |
|
|
226 |
|
|
227 return tree_constant (m); |
|
|
228 } |
|
|
229 |
|
|
230 tree_constant |
|
164
|
231 identity_matrix (const tree_constant& a) |
|
1
|
232 { |
|
|
233 int nr, nc; |
|
143
|
234 get_dimensions (a, "eye", nr, nc); // May set error_state. |
|
|
235 |
|
|
236 if (error_state) |
|
|
237 return tree_constant (); |
|
1
|
238 |
|
|
239 Matrix m (nr, nc, 0.0); |
|
|
240 |
|
|
241 if (nr > 0 && nc > 0) |
|
|
242 { |
|
|
243 int n = MIN (nr, nc); |
|
|
244 for (int i = 0; i < n; i++) |
|
|
245 m.elem (i, i) = 1.0; |
|
|
246 } |
|
|
247 |
|
|
248 return tree_constant (m); |
|
|
249 } |
|
|
250 |
|
|
251 tree_constant |
|
164
|
252 identity_matrix (const tree_constant& a, const tree_constant& b) |
|
1
|
253 { |
|
|
254 int nr, nc; |
|
143
|
255 get_dimensions (a, b, "eye", nr, nc); // May set error_state. |
|
|
256 |
|
|
257 if (error_state) |
|
|
258 return tree_constant (); |
|
1
|
259 |
|
|
260 Matrix m (nr, nc, 0.0); |
|
|
261 |
|
|
262 if (nr > 0 && nc > 0) |
|
|
263 { |
|
|
264 int n = MIN (nr, nc); |
|
|
265 for (int i = 0; i < n; i++) |
|
|
266 m.elem (i, i) = 1.0; |
|
|
267 } |
|
|
268 |
|
|
269 return tree_constant (m); |
|
|
270 } |
|
|
271 |
|
|
272 static tree_constant |
|
|
273 find_nonzero_elem_idx (const Matrix& m) |
|
|
274 { |
|
|
275 int count = 0; |
|
|
276 int m_nr = m.rows (); |
|
|
277 int m_nc = m.columns (); |
|
|
278 |
|
|
279 int i; |
|
|
280 for (int j = 0; j < m_nc; j++) |
|
|
281 for (i = 0; i < m_nr; i++) |
|
|
282 if (m.elem (i, j) != 0) |
|
|
283 count++; |
|
|
284 |
|
|
285 Matrix result; |
|
|
286 |
|
|
287 if (count == 0) |
|
|
288 return result; |
|
|
289 |
|
|
290 if (m_nr == 1) |
|
|
291 { |
|
|
292 result.resize (1, count); |
|
|
293 count = 0; |
|
|
294 for (j = 0; j < m_nc; j++) |
|
|
295 if (m.elem (0, j) != 0) |
|
|
296 { |
|
|
297 result (0, count) = j + 1; |
|
|
298 count++; |
|
|
299 } |
|
|
300 return tree_constant (result); |
|
|
301 } |
|
|
302 else |
|
|
303 { |
|
|
304 ColumnVector v (count); |
|
|
305 count = 0; |
|
|
306 for (j = 0; j < m_nc; j++) |
|
|
307 for (i = 0; i < m_nr; i++) |
|
|
308 if (m.elem (i, j) != 0) |
|
|
309 { |
|
|
310 v.elem (count) = m_nr * j + i + 1; |
|
|
311 count++; |
|
|
312 } |
|
|
313 return tree_constant (v, 1); // Always make a column vector. |
|
|
314 } |
|
|
315 } |
|
|
316 |
|
|
317 static tree_constant |
|
|
318 find_nonzero_elem_idx (const ComplexMatrix& m) |
|
|
319 { |
|
|
320 int count = 0; |
|
|
321 int m_nr = m.rows (); |
|
|
322 int m_nc = m.columns (); |
|
|
323 |
|
|
324 for (int j = 0; j < m_nc; j++) |
|
|
325 { |
|
|
326 for (int i = 0; i < m_nr; i++) |
|
|
327 if (m.elem (i, j) != 0) |
|
|
328 count++; |
|
|
329 } |
|
|
330 |
|
|
331 Matrix result; |
|
|
332 |
|
|
333 if (count == 0) |
|
|
334 return result; |
|
|
335 |
|
|
336 if (m_nr == 1) |
|
|
337 { |
|
|
338 result.resize (1, count); |
|
|
339 count = 0; |
|
|
340 for (j = 0; j < m_nc; j++) |
|
|
341 if (m.elem (0, j) != 0) |
|
|
342 { |
|
|
343 result (0, count) = j + 1; |
|
|
344 count++; |
|
|
345 } |
|
|
346 return tree_constant (result); |
|
|
347 } |
|
|
348 else |
|
|
349 { |
|
|
350 ColumnVector v (count); |
|
|
351 count = 0; |
|
|
352 for (j = 0; j < m_nc; j++) |
|
|
353 { |
|
|
354 for (int i = 0; i < m_nr; i++) |
|
|
355 if (m.elem (i, j) != 0) |
|
|
356 { |
|
|
357 v.elem (count) = m_nr * j + i + 1; |
|
|
358 count++; |
|
|
359 } |
|
|
360 } |
|
|
361 return tree_constant (v, 1); // Always make a column vector. |
|
|
362 } |
|
|
363 } |
|
|
364 |
|
|
365 tree_constant |
|
164
|
366 find_nonzero_elem_idx (const tree_constant& a) |
|
1
|
367 { |
|
|
368 tree_constant retval; |
|
|
369 |
|
|
370 tree_constant tmp = a.make_numeric (); |
|
|
371 |
|
|
372 Matrix result; |
|
|
373 |
|
|
374 switch (tmp.const_type ()) |
|
|
375 { |
|
|
376 case tree_constant_rep::matrix_constant: |
|
|
377 { |
|
|
378 Matrix m = tmp.matrix_value (); |
|
|
379 return find_nonzero_elem_idx (m); |
|
|
380 } |
|
|
381 break; |
|
|
382 case tree_constant_rep::scalar_constant: |
|
|
383 { |
|
|
384 double d = tmp.double_value (); |
|
|
385 if (d != 0.0) |
|
|
386 return tree_constant (1.0); |
|
|
387 else |
|
|
388 return tree_constant (result); |
|
|
389 } |
|
|
390 break; |
|
|
391 case tree_constant_rep::complex_matrix_constant: |
|
|
392 { |
|
|
393 ComplexMatrix m = tmp.complex_matrix_value (); |
|
|
394 return find_nonzero_elem_idx (m); |
|
|
395 } |
|
|
396 break; |
|
|
397 case tree_constant_rep::complex_scalar_constant: |
|
|
398 { |
|
|
399 Complex c = tmp.complex_value (); |
|
|
400 if (c != 0.0) |
|
|
401 return tree_constant (1.0); |
|
|
402 else |
|
|
403 return tree_constant (result); |
|
|
404 } |
|
|
405 break; |
|
|
406 default: |
|
|
407 break; |
|
|
408 } |
|
|
409 return retval; |
|
|
410 } |
|
|
411 |
|
51
|
412 // XXX FIXME XXX -- the next two functions (and expm) should really be just |
|
1
|
413 // one... |
|
|
414 |
|
|
415 tree_constant * |
|
164
|
416 matrix_log (const tree_constant& a) |
|
1
|
417 { |
|
|
418 tree_constant *retval = new tree_constant [2]; |
|
|
419 |
|
|
420 tree_constant tmp = a.make_numeric ();; |
|
|
421 |
|
|
422 if (tmp.rows () == 0 || tmp.columns () == 0) |
|
|
423 { |
|
|
424 int flag = user_pref.propagate_empty_matrices; |
|
|
425 if (flag != 0) |
|
|
426 { |
|
|
427 if (flag < 0) |
|
|
428 gripe_empty_arg ("logm", 0); |
|
|
429 Matrix m; |
|
|
430 retval = new tree_constant [2]; |
|
|
431 retval[0] = tree_constant (m); |
|
|
432 return retval; |
|
|
433 } |
|
|
434 else |
|
|
435 gripe_empty_arg ("logm", 1); |
|
|
436 } |
|
|
437 |
|
|
438 switch (tmp.const_type ()) |
|
|
439 { |
|
|
440 case tree_constant_rep::matrix_constant: |
|
|
441 { |
|
|
442 Matrix m = tmp.matrix_value (); |
|
|
443 |
|
|
444 int nr = m.rows (); |
|
|
445 int nc = m.columns (); |
|
|
446 |
|
|
447 if (nr == 0 || nc == 0 || nr != nc) |
|
|
448 gripe_square_matrix_required ("logm"); |
|
|
449 else |
|
|
450 { |
|
|
451 EIG m_eig (m); |
|
|
452 ComplexColumnVector lambda (m_eig.eigenvalues ()); |
|
|
453 ComplexMatrix Q (m_eig.eigenvectors ()); |
|
|
454 |
|
|
455 for (int i = 0; i < nr; i++) |
|
|
456 { |
|
|
457 Complex elt = lambda.elem (i); |
|
|
458 if (imag (elt) == 0.0 && real (elt) > 0.0) |
|
|
459 lambda.elem (i) = log (real (elt)); |
|
|
460 else |
|
|
461 lambda.elem (i) = log (elt); |
|
|
462 } |
|
|
463 |
|
|
464 ComplexDiagMatrix D (lambda); |
|
|
465 ComplexMatrix result = Q * D * Q.inverse (); |
|
|
466 |
|
|
467 retval[0] = tree_constant (result); |
|
|
468 } |
|
|
469 } |
|
|
470 break; |
|
|
471 case tree_constant_rep::complex_matrix_constant: |
|
|
472 { |
|
|
473 ComplexMatrix m = tmp.complex_matrix_value (); |
|
|
474 |
|
|
475 int nr = m.rows (); |
|
|
476 int nc = m.columns (); |
|
|
477 |
|
|
478 if (nr == 0 || nc == 0 || nr != nc) |
|
|
479 gripe_square_matrix_required ("logm"); |
|
|
480 else |
|
|
481 { |
|
|
482 EIG m_eig (m); |
|
|
483 ComplexColumnVector lambda (m_eig.eigenvalues ()); |
|
|
484 ComplexMatrix Q (m_eig.eigenvectors ()); |
|
|
485 |
|
|
486 for (int i = 0; i < nr; i++) |
|
|
487 { |
|
|
488 Complex elt = lambda.elem (i); |
|
|
489 if (imag (elt) == 0.0 && real (elt) > 0.0) |
|
|
490 lambda.elem (i) = log (real (elt)); |
|
|
491 else |
|
|
492 lambda.elem (i) = log (elt); |
|
|
493 } |
|
|
494 |
|
|
495 ComplexDiagMatrix D (lambda); |
|
|
496 ComplexMatrix result = Q * D * Q.inverse (); |
|
|
497 |
|
|
498 retval[0] = tree_constant (result); |
|
|
499 } |
|
|
500 } |
|
|
501 break; |
|
|
502 case tree_constant_rep::scalar_constant: |
|
|
503 { |
|
|
504 double d = tmp.double_value (); |
|
|
505 if (d > 0.0) |
|
|
506 retval[0] = tree_constant (log (d)); |
|
|
507 else |
|
|
508 { |
|
|
509 Complex dtmp (d); |
|
|
510 retval[0] = tree_constant (log (dtmp)); |
|
|
511 } |
|
|
512 } |
|
|
513 break; |
|
|
514 case tree_constant_rep::complex_scalar_constant: |
|
|
515 { |
|
|
516 Complex c = tmp.complex_value (); |
|
|
517 retval[0] = tree_constant (log (c)); |
|
|
518 } |
|
|
519 break; |
|
|
520 default: |
|
|
521 break; |
|
|
522 } |
|
|
523 return retval; |
|
|
524 } |
|
|
525 |
|
|
526 tree_constant * |
|
164
|
527 matrix_sqrt (const tree_constant& a) |
|
1
|
528 { |
|
|
529 tree_constant *retval = new tree_constant [2]; |
|
|
530 |
|
|
531 tree_constant tmp = a.make_numeric ();; |
|
|
532 |
|
|
533 if (tmp.rows () == 0 || tmp.columns () == 0) |
|
|
534 { |
|
|
535 int flag = user_pref.propagate_empty_matrices; |
|
|
536 if (flag != 0) |
|
|
537 { |
|
|
538 if (flag < 0) |
|
|
539 gripe_empty_arg ("sqrtm", 0); |
|
|
540 Matrix m; |
|
|
541 retval = new tree_constant [2]; |
|
|
542 retval[0] = tree_constant (m); |
|
|
543 return retval; |
|
|
544 } |
|
|
545 else |
|
|
546 gripe_empty_arg ("sqrtm", 1); |
|
|
547 } |
|
|
548 |
|
|
549 switch (tmp.const_type ()) |
|
|
550 { |
|
|
551 case tree_constant_rep::matrix_constant: |
|
|
552 { |
|
|
553 Matrix m = tmp.matrix_value (); |
|
|
554 |
|
|
555 int nr = m.rows (); |
|
|
556 int nc = m.columns (); |
|
|
557 |
|
|
558 if (nr == 0 || nc == 0 || nr != nc) |
|
|
559 gripe_square_matrix_required ("sqrtm"); |
|
|
560 else |
|
|
561 { |
|
|
562 EIG m_eig (m); |
|
|
563 ComplexColumnVector lambda (m_eig.eigenvalues ()); |
|
|
564 ComplexMatrix Q (m_eig.eigenvectors ()); |
|
|
565 |
|
|
566 for (int i = 0; i < nr; i++) |
|
|
567 { |
|
|
568 Complex elt = lambda.elem (i); |
|
|
569 if (imag (elt) == 0.0 && real (elt) > 0.0) |
|
|
570 lambda.elem (i) = sqrt (real (elt)); |
|
|
571 else |
|
|
572 lambda.elem (i) = sqrt (elt); |
|
|
573 } |
|
|
574 |
|
|
575 ComplexDiagMatrix D (lambda); |
|
|
576 ComplexMatrix result = Q * D * Q.inverse (); |
|
|
577 |
|
|
578 retval[0] = tree_constant (result); |
|
|
579 } |
|
|
580 } |
|
|
581 break; |
|
|
582 case tree_constant_rep::complex_matrix_constant: |
|
|
583 { |
|
|
584 ComplexMatrix m = tmp.complex_matrix_value (); |
|
|
585 |
|
|
586 int nr = m.rows (); |
|
|
587 int nc = m.columns (); |
|
|
588 |
|
|
589 if (nr == 0 || nc == 0 || nr != nc) |
|
|
590 gripe_square_matrix_required ("sqrtm"); |
|
|
591 else |
|
|
592 { |
|
|
593 EIG m_eig (m); |
|
|
594 ComplexColumnVector lambda (m_eig.eigenvalues ()); |
|
|
595 ComplexMatrix Q (m_eig.eigenvectors ()); |
|
|
596 |
|
|
597 for (int i = 0; i < nr; i++) |
|
|
598 { |
|
|
599 Complex elt = lambda.elem (i); |
|
|
600 if (imag (elt) == 0.0 && real (elt) > 0.0) |
|
|
601 lambda.elem (i) = sqrt (real (elt)); |
|
|
602 else |
|
|
603 lambda.elem (i) = sqrt (elt); |
|
|
604 } |
|
|
605 |
|
|
606 ComplexDiagMatrix D (lambda); |
|
|
607 ComplexMatrix result = Q * D * Q.inverse (); |
|
|
608 |
|
|
609 retval[0] = tree_constant (result); |
|
|
610 } |
|
|
611 } |
|
|
612 break; |
|
|
613 case tree_constant_rep::scalar_constant: |
|
|
614 { |
|
|
615 double d = tmp.double_value (); |
|
|
616 if (d > 0.0) |
|
|
617 retval[0] = tree_constant (sqrt (d)); |
|
|
618 else |
|
|
619 { |
|
|
620 Complex dtmp (d); |
|
|
621 retval[0] = tree_constant (sqrt (dtmp)); |
|
|
622 } |
|
|
623 } |
|
|
624 break; |
|
|
625 case tree_constant_rep::complex_scalar_constant: |
|
|
626 { |
|
|
627 Complex c = tmp.complex_value (); |
|
|
628 retval[0] = tree_constant (log (c)); |
|
|
629 } |
|
|
630 break; |
|
|
631 default: |
|
|
632 break; |
|
|
633 } |
|
|
634 return retval; |
|
|
635 } |
|
|
636 |
|
|
637 tree_constant * |
|
164
|
638 column_max (const tree_constant *args, int nargin, int nargout) |
|
1
|
639 { |
|
|
640 tree_constant *retval = NULL_TREE_CONST; |
|
|
641 |
|
|
642 tree_constant arg1; |
|
|
643 tree_constant arg2; |
|
|
644 tree_constant_rep::constant_type arg1_type = |
|
|
645 tree_constant_rep::unknown_constant; |
|
|
646 tree_constant_rep::constant_type arg2_type = |
|
|
647 tree_constant_rep::unknown_constant; |
|
|
648 |
|
|
649 switch (nargin) |
|
|
650 { |
|
|
651 case 3: |
|
|
652 arg2 = args[2].make_numeric (); |
|
|
653 arg2_type = arg2.const_type (); |
|
|
654 // Fall through... |
|
|
655 case 2: |
|
|
656 arg1 = args[1].make_numeric (); |
|
|
657 arg1_type = arg1.const_type (); |
|
|
658 break; |
|
|
659 default: |
|
|
660 panic_impossible (); |
|
|
661 break; |
|
|
662 } |
|
|
663 |
|
|
664 if (nargin == 2 && nargout == 1) |
|
|
665 { |
|
|
666 retval = new tree_constant [2]; |
|
199
|
667 switch (arg1_type) |
|
1
|
668 { |
|
199
|
669 case tree_constant_rep::scalar_constant: |
|
|
670 retval[0] = tree_constant (arg1.double_value ()); |
|
|
671 break; |
|
|
672 case tree_constant_rep::complex_scalar_constant: |
|
|
673 retval[0] = tree_constant (arg1.complex_value ()); |
|
|
674 break; |
|
|
675 case tree_constant_rep::matrix_constant: |
|
|
676 { |
|
|
677 Matrix m = arg1.matrix_value (); |
|
|
678 if (m.rows () == 1) |
|
|
679 retval[0] = tree_constant (m.row_max ()); |
|
|
680 else |
|
|
681 retval[0] = tree_constant (m.column_max (), 0); |
|
|
682 } |
|
|
683 break; |
|
|
684 case tree_constant_rep::complex_matrix_constant: |
|
|
685 { |
|
|
686 ComplexMatrix m = arg1.complex_matrix_value (); |
|
|
687 if (m.rows () == 1) |
|
|
688 retval[0] = tree_constant (m.row_max ()); |
|
|
689 else |
|
|
690 retval[0] = tree_constant (m.column_max (), 0); |
|
|
691 } |
|
|
692 break; |
|
|
693 default: |
|
|
694 panic_impossible (); |
|
|
695 break; |
|
1
|
696 } |
|
|
697 } |
|
|
698 else if (nargin == 2 && nargout == 2) |
|
199
|
699 { |
|
|
700 retval = new tree_constant [2]; |
|
|
701 switch (arg1_type) |
|
|
702 { |
|
|
703 case tree_constant_rep::scalar_constant: |
|
|
704 { |
|
|
705 retval[0] = tree_constant (arg1.double_value ()); |
|
|
706 retval[1] = tree_constant (1); |
|
|
707 } |
|
|
708 break; |
|
|
709 case tree_constant_rep::complex_scalar_constant: |
|
|
710 { |
|
|
711 retval[0] = tree_constant (arg1.complex_value ()); |
|
|
712 retval[1] = tree_constant (1); |
|
|
713 } |
|
|
714 break; |
|
|
715 case tree_constant_rep::matrix_constant: |
|
|
716 { |
|
|
717 Matrix m = arg1.matrix_value (); |
|
|
718 if (m.rows () == 1) |
|
|
719 { |
|
|
720 retval[0] = tree_constant (m.row_max ()); |
|
|
721 retval[1] = tree_constant (m.row_max_loc ()); |
|
|
722 } |
|
|
723 else |
|
|
724 { |
|
|
725 retval[0] = tree_constant (m.column_max (), 0); |
|
|
726 retval[1] = tree_constant (m.column_max_loc (), 0); |
|
|
727 } |
|
|
728 } |
|
|
729 break; |
|
|
730 case tree_constant_rep::complex_matrix_constant: |
|
|
731 { |
|
|
732 ComplexMatrix m = arg1.complex_matrix_value (); |
|
|
733 if (m.rows () == 1) |
|
|
734 { |
|
|
735 retval[0] = tree_constant (m.row_max ()); |
|
|
736 retval[1] = tree_constant (m.row_max_loc ()); |
|
|
737 } |
|
|
738 else |
|
|
739 { |
|
|
740 retval[0] = tree_constant (m.column_max (), 0); |
|
|
741 retval[1] = tree_constant (m.column_max_loc (), 0); |
|
|
742 } |
|
|
743 } |
|
|
744 break; |
|
|
745 default: |
|
|
746 panic_impossible (); |
|
|
747 break; |
|
|
748 } |
|
|
749 } |
|
1
|
750 else if (nargin == 3) |
|
|
751 { |
|
|
752 if (arg1.rows () == arg2.rows () |
|
|
753 && arg1.columns () == arg2.columns ()) |
|
|
754 { |
|
|
755 retval = new tree_constant [2]; |
|
199
|
756 switch (arg1_type) |
|
|
757 { |
|
|
758 case tree_constant_rep::scalar_constant: |
|
|
759 { |
|
|
760 double result; |
|
|
761 double a_elem = arg1.double_value (); |
|
|
762 double b_elem = arg2.double_value (); |
|
433
|
763 result = MAX (a_elem, b_elem); |
|
199
|
764 retval[0] = tree_constant (result); |
|
|
765 } |
|
|
766 break; |
|
|
767 case tree_constant_rep::complex_scalar_constant: |
|
|
768 { |
|
|
769 Complex result; |
|
|
770 Complex a_elem = arg1.complex_value (); |
|
|
771 Complex b_elem = arg2.complex_value (); |
|
433
|
772 if (abs (a_elem) > abs (b_elem)) |
|
199
|
773 result = a_elem; |
|
|
774 else |
|
|
775 result = b_elem; |
|
|
776 retval[0] = tree_constant (result); |
|
|
777 } |
|
|
778 break; |
|
|
779 case tree_constant_rep::matrix_constant: |
|
|
780 { |
|
|
781 Matrix result; |
|
|
782 result = max (arg1.matrix_value (), arg2.matrix_value ()); |
|
|
783 retval[0] = tree_constant (result); |
|
|
784 } |
|
|
785 break; |
|
|
786 case tree_constant_rep::complex_matrix_constant: |
|
|
787 { |
|
|
788 ComplexMatrix result; |
|
|
789 result = max (arg1.complex_matrix_value (), |
|
|
790 arg2.complex_matrix_value ()); |
|
|
791 retval[0] = tree_constant (result); |
|
|
792 } |
|
|
793 break; |
|
|
794 default: |
|
|
795 panic_impossible (); |
|
|
796 break; |
|
1
|
797 } |
|
|
798 } |
|
|
799 else |
|
217
|
800 error ("max: nonconformant matrices"); |
|
1
|
801 } |
|
|
802 else |
|
|
803 panic_impossible (); |
|
|
804 |
|
|
805 return retval; |
|
|
806 } |
|
|
807 |
|
|
808 tree_constant * |
|
164
|
809 column_min (const tree_constant *args, int nargin, int nargout) |
|
1
|
810 { |
|
|
811 tree_constant *retval = NULL_TREE_CONST; |
|
|
812 |
|
|
813 tree_constant arg1; |
|
|
814 tree_constant arg2; |
|
|
815 tree_constant_rep::constant_type arg1_type = |
|
|
816 tree_constant_rep::unknown_constant; |
|
|
817 tree_constant_rep::constant_type arg2_type = |
|
|
818 tree_constant_rep::unknown_constant; |
|
|
819 |
|
|
820 switch (nargin) |
|
|
821 { |
|
|
822 case 3: |
|
|
823 arg2 = args[2].make_numeric (); |
|
|
824 arg2_type = arg2.const_type (); |
|
|
825 // Fall through... |
|
|
826 case 2: |
|
|
827 arg1 = args[1].make_numeric (); |
|
|
828 arg1_type = arg1.const_type (); |
|
|
829 break; |
|
|
830 default: |
|
|
831 panic_impossible (); |
|
|
832 break; |
|
|
833 } |
|
|
834 |
|
|
835 if (nargin == 2 && nargout == 1) |
|
|
836 { |
|
|
837 retval = new tree_constant [2]; |
|
199
|
838 switch (arg1_type) |
|
1
|
839 { |
|
199
|
840 case tree_constant_rep::scalar_constant: |
|
|
841 retval[0] = tree_constant (arg1.double_value ()); |
|
|
842 break; |
|
|
843 case tree_constant_rep::complex_scalar_constant: |
|
|
844 retval[0] = tree_constant (arg1.complex_value ()); |
|
|
845 break; |
|
|
846 case tree_constant_rep::matrix_constant: |
|
|
847 { |
|
|
848 Matrix m = arg1.matrix_value (); |
|
|
849 if (m.rows () == 1) |
|
|
850 retval[0] = tree_constant (m.row_min ()); |
|
|
851 else |
|
|
852 retval[0] = tree_constant (m.column_min (), 0); |
|
|
853 } |
|
|
854 break; |
|
|
855 case tree_constant_rep::complex_matrix_constant: |
|
|
856 { |
|
|
857 ComplexMatrix m = arg1.complex_matrix_value (); |
|
|
858 if (m.rows () == 1) |
|
|
859 retval[0] = tree_constant (m.row_min ()); |
|
|
860 else |
|
|
861 retval[0] = tree_constant (m.column_min (), 0); |
|
|
862 } |
|
|
863 break; |
|
|
864 default: |
|
|
865 panic_impossible (); |
|
|
866 break; |
|
1
|
867 } |
|
|
868 } |
|
|
869 else if (nargin == 2 && nargout == 2) |
|
199
|
870 { |
|
|
871 retval = new tree_constant [2]; |
|
|
872 switch (arg1_type) |
|
|
873 { |
|
|
874 case tree_constant_rep::scalar_constant: |
|
|
875 { |
|
|
876 retval[0] = tree_constant (arg1.double_value ()); |
|
|
877 retval[1] = tree_constant (1); |
|
|
878 } |
|
|
879 break; |
|
|
880 case tree_constant_rep::complex_scalar_constant: |
|
|
881 { |
|
|
882 retval[0] = tree_constant (arg1.complex_value ()); |
|
|
883 retval[1] = tree_constant (1); |
|
|
884 } |
|
|
885 break; |
|
|
886 case tree_constant_rep::matrix_constant: |
|
|
887 { |
|
|
888 Matrix m = arg1.matrix_value (); |
|
|
889 if (m.rows () == 1) |
|
|
890 { |
|
|
891 retval[0] = tree_constant (m.row_min ()); |
|
|
892 retval[1] = tree_constant (m.row_min_loc ()); |
|
|
893 } |
|
|
894 else |
|
|
895 { |
|
|
896 retval[0] = tree_constant (m.column_min (), 0); |
|
|
897 retval[1] = tree_constant (m.column_min_loc (), 0); |
|
|
898 } |
|
|
899 } |
|
|
900 break; |
|
|
901 case tree_constant_rep::complex_matrix_constant: |
|
|
902 { |
|
|
903 ComplexMatrix m = arg1.complex_matrix_value (); |
|
|
904 if (m.rows () == 1) |
|
|
905 { |
|
|
906 retval[0] = tree_constant (m.row_min ()); |
|
|
907 retval[1] = tree_constant (m.row_min_loc ()); |
|
|
908 } |
|
|
909 else |
|
|
910 { |
|
|
911 retval[0] = tree_constant (m.column_min (), 0); |
|
|
912 retval[1] = tree_constant (m.column_min_loc (), 0); |
|
|
913 } |
|
|
914 } |
|
|
915 break; |
|
|
916 default: |
|
|
917 panic_impossible (); |
|
|
918 break; |
|
|
919 } |
|
|
920 } |
|
1
|
921 else if (nargin == 3) |
|
|
922 { |
|
|
923 if (arg1.rows () == arg2.rows () |
|
|
924 && arg1.columns () == arg2.columns ()) |
|
|
925 { |
|
|
926 retval = new tree_constant [2]; |
|
199
|
927 switch (arg1_type) |
|
|
928 { |
|
|
929 case tree_constant_rep::scalar_constant: |
|
|
930 { |
|
|
931 double result; |
|
|
932 double a_elem = arg1.double_value (); |
|
|
933 double b_elem = arg2.double_value (); |
|
|
934 result = MIN (a_elem, b_elem); |
|
|
935 retval[0] = tree_constant (result); |
|
|
936 } |
|
|
937 break; |
|
|
938 case tree_constant_rep::complex_scalar_constant: |
|
|
939 { |
|
|
940 Complex result; |
|
|
941 Complex a_elem = arg1.complex_value (); |
|
|
942 Complex b_elem = arg2.complex_value (); |
|
433
|
943 if (abs (a_elem) < abs (b_elem)) |
|
199
|
944 result = a_elem; |
|
|
945 else |
|
|
946 result = b_elem; |
|
|
947 retval[0] = tree_constant (result); |
|
|
948 } |
|
|
949 break; |
|
|
950 case tree_constant_rep::matrix_constant: |
|
|
951 { |
|
|
952 Matrix result; |
|
|
953 result = min (arg1.matrix_value (), arg2.matrix_value ()); |
|
|
954 retval[0] = tree_constant (result); |
|
|
955 } |
|
|
956 break; |
|
|
957 case tree_constant_rep::complex_matrix_constant: |
|
|
958 { |
|
|
959 ComplexMatrix result; |
|
|
960 result = min (arg1.complex_matrix_value (), |
|
|
961 arg2.complex_matrix_value ()); |
|
|
962 retval[0] = tree_constant (result); |
|
|
963 } |
|
|
964 break; |
|
|
965 default: |
|
|
966 panic_impossible (); |
|
|
967 break; |
|
1
|
968 } |
|
|
969 } |
|
|
970 else |
|
217
|
971 error ("min: nonconformant matrices"); |
|
1
|
972 } |
|
|
973 else |
|
|
974 panic_impossible (); |
|
|
975 |
|
|
976 return retval; |
|
|
977 } |
|
|
978 |
|
|
979 static void |
|
|
980 mx_sort (Matrix& m, Matrix& idx, int return_idx) |
|
|
981 { |
|
|
982 int nr = m.rows (); |
|
|
983 int nc = m.columns (); |
|
|
984 idx.resize (nr, nc); |
|
|
985 int i, j; |
|
|
986 |
|
|
987 if (return_idx) |
|
|
988 { |
|
|
989 for (j = 0; j < nc; j++) |
|
|
990 for (i = 0; i < nr; i++) |
|
|
991 idx.elem (i, j) = i+1; |
|
|
992 } |
|
|
993 |
|
|
994 for (j = 0; j < nc; j++) |
|
|
995 { |
|
|
996 for (int gap = nr/2; gap > 0; gap /= 2) |
|
|
997 for (i = gap; i < nr; i++) |
|
|
998 for (int k = i - gap; |
|
|
999 k >= 0 && m.elem (k, j) > m.elem (k+gap, j); |
|
|
1000 k -= gap) |
|
|
1001 { |
|
|
1002 double tmp = m.elem (k, j); |
|
|
1003 m.elem (k, j) = m.elem (k+gap, j); |
|
|
1004 m.elem (k+gap, j) = tmp; |
|
|
1005 |
|
|
1006 if (return_idx) |
|
|
1007 { |
|
|
1008 double tmp = idx.elem (k, j); |
|
|
1009 idx.elem (k, j) = idx.elem (k+gap, j); |
|
|
1010 idx.elem (k+gap, j) = tmp; |
|
|
1011 } |
|
|
1012 } |
|
|
1013 } |
|
|
1014 } |
|
|
1015 |
|
|
1016 static void |
|
|
1017 mx_sort (RowVector& v, RowVector& idx, int return_idx) |
|
|
1018 { |
|
|
1019 int n = v.capacity (); |
|
|
1020 idx.resize (n); |
|
|
1021 int i; |
|
|
1022 |
|
|
1023 if (return_idx) |
|
|
1024 for (i = 0; i < n; i++) |
|
|
1025 idx.elem (i) = i+1; |
|
|
1026 |
|
|
1027 for (int gap = n/2; gap > 0; gap /= 2) |
|
|
1028 for (i = gap; i < n; i++) |
|
|
1029 for (int k = i - gap; |
|
|
1030 k >= 0 && v.elem (k) > v.elem (k+gap); |
|
|
1031 k -= gap) |
|
|
1032 { |
|
|
1033 double tmp = v.elem (k); |
|
|
1034 v.elem (k) = v.elem (k+gap); |
|
|
1035 v.elem (k+gap) = tmp; |
|
|
1036 |
|
|
1037 if (return_idx) |
|
|
1038 { |
|
|
1039 double tmp = idx.elem (k); |
|
|
1040 idx.elem (k) = idx.elem (k+gap); |
|
|
1041 idx.elem (k+gap) = tmp; |
|
|
1042 } |
|
|
1043 } |
|
|
1044 } |
|
|
1045 |
|
|
1046 static void |
|
|
1047 mx_sort (ComplexMatrix& cm, Matrix& idx, int return_idx) |
|
|
1048 { |
|
|
1049 int nr = cm.rows (); |
|
|
1050 int nc = cm.columns (); |
|
|
1051 idx.resize (nr, nc); |
|
|
1052 int i, j; |
|
|
1053 |
|
|
1054 if (return_idx) |
|
|
1055 { |
|
|
1056 for (j = 0; j < nc; j++) |
|
|
1057 for (i = 0; i < nr; i++) |
|
|
1058 idx.elem (i, j) = i+1; |
|
|
1059 } |
|
|
1060 |
|
|
1061 for (j = 0; j < nc; j++) |
|
|
1062 { |
|
|
1063 for (int gap = nr/2; gap > 0; gap /= 2) |
|
|
1064 for (i = gap; i < nr; i++) |
|
|
1065 for (int k = i - gap; |
|
|
1066 k >= 0 && abs (cm.elem (k, j)) > abs (cm.elem (k+gap, j)); |
|
|
1067 k -= gap) |
|
|
1068 { |
|
|
1069 Complex ctmp = cm.elem (k, j); |
|
|
1070 cm.elem (k, j) = cm.elem (k+gap, j); |
|
|
1071 cm.elem (k+gap, j) = ctmp; |
|
|
1072 |
|
|
1073 if (return_idx) |
|
|
1074 { |
|
|
1075 double tmp = idx.elem (k, j); |
|
|
1076 idx.elem (k, j) = idx.elem (k+gap, j); |
|
|
1077 idx.elem (k+gap, j) = tmp; |
|
|
1078 } |
|
|
1079 } |
|
|
1080 } |
|
|
1081 } |
|
|
1082 |
|
|
1083 static void |
|
|
1084 mx_sort (ComplexRowVector& cv, RowVector& idx, int return_idx) |
|
|
1085 { |
|
|
1086 int n = cv.capacity (); |
|
|
1087 idx.resize (n); |
|
|
1088 int i; |
|
|
1089 |
|
|
1090 if (return_idx) |
|
|
1091 for (i = 0; i < n; i++) |
|
|
1092 idx.elem (i) = i+1; |
|
|
1093 |
|
|
1094 for (int gap = n/2; gap > 0; gap /= 2) |
|
|
1095 for (i = gap; i < n; i++) |
|
|
1096 for (int k = i - gap; |
|
|
1097 k >= 0 && abs (cv.elem (k)) > abs (cv.elem (k+gap)); |
|
|
1098 k -= gap) |
|
|
1099 { |
|
|
1100 Complex tmp = cv.elem (k); |
|
|
1101 cv.elem (k) = cv.elem (k+gap); |
|
|
1102 cv.elem (k+gap) = tmp; |
|
|
1103 |
|
|
1104 if (return_idx) |
|
|
1105 { |
|
|
1106 double tmp = idx.elem (k); |
|
|
1107 idx.elem (k) = idx.elem (k+gap); |
|
|
1108 idx.elem (k+gap) = tmp; |
|
|
1109 } |
|
|
1110 } |
|
|
1111 } |
|
|
1112 |
|
|
1113 tree_constant * |
|
164
|
1114 sort (const tree_constant *args, int nargin, int nargout) |
|
1
|
1115 { |
|
|
1116 // Assumes that we have been given the correct number of arguments. |
|
|
1117 |
|
|
1118 tree_constant *retval = NULL_TREE_CONST; |
|
|
1119 |
|
|
1120 int return_idx = nargout > 1; |
|
|
1121 if (return_idx) |
|
|
1122 retval = new tree_constant [3]; |
|
|
1123 else |
|
|
1124 retval = new tree_constant [2]; |
|
|
1125 |
|
|
1126 switch (args[1].const_type ()) |
|
|
1127 { |
|
|
1128 case tree_constant_rep::scalar_constant: |
|
|
1129 { |
|
|
1130 retval [0] = tree_constant (args[1].double_value ()); |
|
|
1131 if (return_idx) |
|
|
1132 retval [1] = tree_constant (1.0); |
|
|
1133 } |
|
|
1134 break; |
|
|
1135 case tree_constant_rep::complex_scalar_constant: |
|
|
1136 { |
|
|
1137 retval [0] = tree_constant (args[1].complex_value ()); |
|
|
1138 if (return_idx) |
|
|
1139 retval [1] = tree_constant (1.0); |
|
|
1140 } |
|
|
1141 break; |
|
|
1142 case tree_constant_rep::string_constant: |
|
|
1143 case tree_constant_rep::range_constant: |
|
|
1144 case tree_constant_rep::matrix_constant: |
|
|
1145 { |
|
|
1146 Matrix m = args[1].to_matrix (); |
|
|
1147 if (m.rows () == 1) |
|
|
1148 { |
|
|
1149 int nc = m.columns (); |
|
|
1150 RowVector v (nc); |
|
|
1151 for (int i = 0; i < nc; i++) |
|
|
1152 v.elem (i) = m.elem (0, i); |
|
|
1153 RowVector idx; |
|
|
1154 mx_sort (v, idx, return_idx); |
|
|
1155 |
|
|
1156 retval [0] = tree_constant (v, 0); |
|
|
1157 if (return_idx) |
|
|
1158 retval [1] = tree_constant (idx, 0); |
|
|
1159 } |
|
|
1160 else |
|
|
1161 { |
|
|
1162 // Sorts m in place, optionally computes index Matrix. |
|
|
1163 Matrix idx; |
|
|
1164 mx_sort (m, idx, return_idx); |
|
|
1165 |
|
|
1166 retval [0] = tree_constant (m); |
|
|
1167 if (return_idx) |
|
|
1168 retval [1] = tree_constant (idx); |
|
|
1169 } |
|
|
1170 } |
|
|
1171 break; |
|
|
1172 case tree_constant_rep::complex_matrix_constant: |
|
|
1173 { |
|
|
1174 ComplexMatrix cm = args[1].complex_matrix_value (); |
|
|
1175 if (cm.rows () == 1) |
|
|
1176 { |
|
|
1177 int nc = cm.columns (); |
|
|
1178 ComplexRowVector cv (nc); |
|
|
1179 for (int i = 0; i < nc; i++) |
|
|
1180 cv.elem (i) = cm.elem (0, i); |
|
|
1181 RowVector idx; |
|
|
1182 mx_sort (cv, idx, return_idx); |
|
|
1183 |
|
|
1184 retval [0] = tree_constant (cv, 0); |
|
|
1185 if (return_idx) |
|
|
1186 retval [1] = tree_constant (idx, 0); |
|
|
1187 } |
|
|
1188 else |
|
|
1189 { |
|
|
1190 // Sorts cm in place, optionally computes index Matrix. |
|
|
1191 Matrix idx; |
|
|
1192 mx_sort (cm, idx, return_idx); |
|
|
1193 |
|
|
1194 retval [0] = tree_constant (cm); |
|
|
1195 if (return_idx) |
|
|
1196 retval [1] = tree_constant (idx); |
|
|
1197 } |
|
|
1198 } |
|
|
1199 break; |
|
|
1200 default: |
|
|
1201 panic_impossible (); |
|
|
1202 break; |
|
|
1203 } |
|
|
1204 |
|
|
1205 return retval; |
|
|
1206 } |
|
|
1207 |
|
|
1208 tree_constant * |
|
164
|
1209 feval (const tree_constant *args, int nargin, int nargout) |
|
1
|
1210 { |
|
|
1211 // Assumes that we have been given the correct number of arguments. |
|
|
1212 |
|
|
1213 tree_constant *retval = NULL_TREE_CONST; |
|
|
1214 |
|
|
1215 tree *fcn = is_valid_function (args[1], "feval", 1); |
|
|
1216 if (fcn != NULL_TREE) |
|
|
1217 { |
|
|
1218 args++; |
|
|
1219 nargin--; |
|
|
1220 if (nargin > 1) |
|
|
1221 retval = fcn->eval (args, nargin, nargout, 0); |
|
|
1222 else |
|
|
1223 retval = fcn->eval (0, nargout); |
|
|
1224 } |
|
|
1225 |
|
|
1226 return retval; |
|
|
1227 } |
|
|
1228 |
|
|
1229 tree_constant |
|
164
|
1230 eval_string (const char *string, int print, int ans_assign, |
|
1
|
1231 int& parse_status) |
|
|
1232 { |
|
|
1233 begin_unwind_frame ("eval_string"); |
|
|
1234 |
|
|
1235 unwind_protect_int (get_input_from_eval_string); |
|
|
1236 unwind_protect_ptr (global_command); |
|
|
1237 unwind_protect_ptr (current_eval_string); |
|
|
1238 |
|
|
1239 get_input_from_eval_string = 1; |
|
|
1240 current_eval_string = string; |
|
|
1241 |
|
|
1242 YY_BUFFER_STATE old_buf = current_buffer (); |
|
|
1243 YY_BUFFER_STATE new_buf = create_buffer ((FILE *) NULL); |
|
|
1244 |
|
|
1245 add_unwind_protect (restore_input_buffer, (void *) old_buf); |
|
|
1246 add_unwind_protect (delete_input_buffer, (void *) new_buf); |
|
|
1247 |
|
|
1248 switch_to_buffer (new_buf); |
|
|
1249 |
|
|
1250 unwind_protect_ptr (curr_sym_tab); |
|
287
|
1251 |
|
|
1252 reset_parser (); |
|
1
|
1253 |
|
|
1254 parse_status = yyparse (); |
|
|
1255 |
|
|
1256 // Important to reset the idea of where input is coming from before |
|
|
1257 // trying to eval the command we just parsed -- it might contain the |
|
338
|
1258 // name of an function file that still needs to be parsed! |
|
1
|
1259 |
|
|
1260 tree *command = global_command; |
|
|
1261 |
|
|
1262 run_unwind_frame ("eval_string"); |
|
|
1263 |
|
|
1264 tree_constant retval; |
|
|
1265 |
|
|
1266 if (parse_status == 0 && command != NULL_TREE) |
|
|
1267 { |
|
|
1268 retval = command->eval (print); |
|
|
1269 delete command; |
|
|
1270 } |
|
|
1271 |
|
|
1272 return retval; |
|
|
1273 } |
|
|
1274 |
|
|
1275 tree_constant |
|
164
|
1276 eval_string (const tree_constant& arg, int& parse_status) |
|
1
|
1277 { |
|
|
1278 if (! arg.is_string_type ()) |
|
|
1279 { |
|
|
1280 error ("eval: expecting string argument"); |
|
|
1281 return -1; |
|
|
1282 } |
|
|
1283 |
|
|
1284 char *string = arg.string_value (); |
|
|
1285 |
|
|
1286 // Yes Virginia, we always print here... |
|
|
1287 |
|
|
1288 return eval_string (string, 1, 1, parse_status); |
|
|
1289 } |
|
|
1290 |
|
|
1291 static int |
|
164
|
1292 match_sans_spaces (const char *standard, const char *test) |
|
1
|
1293 { |
|
164
|
1294 const char *tp = test; |
|
1
|
1295 while (*tp == ' ' || *tp == '\t') |
|
|
1296 tp++; |
|
|
1297 |
|
164
|
1298 const char *ep = test + strlen (test) - 1; |
|
1
|
1299 while (*ep == ' ' || *ep == '\t') |
|
|
1300 ep--; |
|
|
1301 |
|
|
1302 int len = ep - tp + 1; |
|
|
1303 |
|
|
1304 return (strncmp (standard, tp, len) == 0); |
|
|
1305 } |
|
|
1306 |
|
|
1307 tree_constant |
|
164
|
1308 get_user_input (const tree_constant *args, int nargin, int nargout, |
|
|
1309 int debug = 0) |
|
1
|
1310 { |
|
|
1311 tree_constant retval; |
|
|
1312 |
|
|
1313 int read_as_string = 0; |
|
|
1314 if (nargin == 3) |
|
|
1315 { |
|
|
1316 if (args[2].is_string_type () |
|
|
1317 && strcmp ("s", args[2].string_value ()) == 0) |
|
|
1318 read_as_string++; |
|
|
1319 else |
|
|
1320 { |
|
|
1321 error ("input: unrecognized second argument"); |
|
|
1322 return retval; |
|
|
1323 } |
|
|
1324 } |
|
|
1325 |
|
|
1326 char *prompt = "debug> "; |
|
|
1327 if (nargin > 1) |
|
|
1328 { |
|
|
1329 if (args[1].is_string_type ()) |
|
|
1330 prompt = args[1].string_value (); |
|
|
1331 else |
|
|
1332 { |
|
|
1333 error ("input: unrecognized argument"); |
|
|
1334 return retval; |
|
|
1335 } |
|
|
1336 } |
|
|
1337 |
|
|
1338 again: |
|
|
1339 |
|
|
1340 flush_output_to_pager (); |
|
|
1341 |
|
|
1342 char *input_buf = gnu_readline (prompt); |
|
|
1343 |
|
|
1344 if (input_buf != (char *) NULL) |
|
|
1345 { |
|
|
1346 if (input_buf) |
|
|
1347 maybe_save_history (input_buf); |
|
|
1348 |
|
|
1349 int len = strlen (input_buf); |
|
|
1350 |
|
|
1351 if (len < 1) |
|
|
1352 { |
|
|
1353 if (debug) |
|
|
1354 goto again; |
|
|
1355 else |
|
|
1356 return retval; |
|
|
1357 } |
|
|
1358 |
|
|
1359 if (match_sans_spaces ("exit", input_buf) |
|
|
1360 || match_sans_spaces ("quit", input_buf) |
|
|
1361 || match_sans_spaces ("return", input_buf)) |
|
|
1362 return tree_constant (); |
|
|
1363 else if (read_as_string) |
|
|
1364 retval = tree_constant (input_buf); |
|
|
1365 else |
|
|
1366 { |
|
|
1367 int parse_status; |
|
|
1368 retval = eval_string (input_buf, 0, 0, parse_status); |
|
|
1369 if (debug && retval.is_defined ()) |
|
|
1370 retval.eval (1); |
|
|
1371 } |
|
|
1372 } |
|
|
1373 else |
|
|
1374 error ("input: reading user-input failed!"); |
|
|
1375 |
|
|
1376 if (debug) |
|
|
1377 goto again; |
|
|
1378 |
|
|
1379 return retval; |
|
|
1380 } |
|
|
1381 |
|
|
1382 /* |
|
|
1383 ;;; Local Variables: *** |
|
|
1384 ;;; mode: C++ *** |
|
|
1385 ;;; page-delimiter: "^/\\*" *** |
|
|
1386 ;;; End: *** |
|
|
1387 */ |
