492
|
1 // The constants for the tree class. -*- C++ -*- |
|
2 /* |
|
3 |
|
4 Copyright (C) 1992, 1993, 1994 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 HAVE_CONFIG_H |
|
25 #include "config.h" |
|
26 #endif |
|
27 |
|
28 #if defined (__GNUG__) |
|
29 #pragma implementation |
|
30 #endif |
|
31 |
|
32 #include <ctype.h> |
|
33 #include <string.h> |
|
34 #include <fstream.h> |
|
35 #include <iostream.h> |
|
36 #include <strstream.h> |
|
37 |
|
38 #include "mx-base.h" |
|
39 #include "Range.h" |
|
40 |
|
41 #include "variables.h" |
|
42 #include "error.h" |
|
43 #include "gripes.h" |
|
44 #include "user-prefs.h" |
|
45 #include "utils.h" |
|
46 #include "pager.h" |
|
47 #include "pr-output.h" |
|
48 #include "tree-const.h" |
|
49 #include "idx-vector.h" |
|
50 |
|
51 #include "tc-inlines.cc" |
|
52 |
|
53 /* |
|
54 * How about a few macros? |
|
55 */ |
|
56 |
|
57 #ifndef MAX |
|
58 #define MAX(a,b) ((a) > (b) ? (a) : (b)) |
|
59 #endif |
|
60 |
|
61 #ifndef MIN |
|
62 #define MIN(a,b) ((a) < (b) ? (a) : (b)) |
|
63 #endif |
|
64 |
|
65 #ifndef ABS |
|
66 #define ABS(x) (((x) < 0) ? (-x) : (x)) |
|
67 #endif |
|
68 |
|
69 /* |
|
70 * The following are used by some of the functions in the |
|
71 * tree_constant_rep class that must deal with real and complex |
|
72 * matrices. This was not done with overloaded or virtual functions |
|
73 * from the Matrix class because there is no clean way to do that -- |
|
74 * the necessary functions (like elem) need to return values of |
|
75 * different types... |
|
76 */ |
|
77 |
|
78 // Given a tree_constant, and the names to be used for the real and |
|
79 // complex matrix and their dimensions, declare a real or complex |
|
80 // matrix, and initialize it from the tree_constant. Note that m, cm, |
|
81 // nr, and nc must not be previously declared, and they must not be |
|
82 // expressions. Since only one of the matrices will be defined after |
|
83 // this macro is used, only one set of dimesions is declared. |
|
84 |
|
85 // This macro only makes sense inside a friend or member function of |
|
86 // the tree_constant_rep class |
|
87 |
|
88 #define REP_RHS_MATRIX(tc,m,cm,nr,nc) \ |
|
89 int nr = 0; \ |
|
90 int nc = 0; \ |
|
91 Matrix m; \ |
|
92 ComplexMatrix cm; \ |
|
93 if ((tc).const_type () == tree_constant_rep::complex_matrix_constant) \ |
|
94 { \ |
|
95 cm = (tc).complex_matrix_value (); \ |
|
96 nr = (cm).rows (); \ |
|
97 nc = (cm).columns (); \ |
|
98 } \ |
|
99 else if ((tc).const_type () == tree_constant_rep::matrix_constant) \ |
|
100 { \ |
|
101 m = (tc).matrix_value (); \ |
|
102 nr = (m).rows (); \ |
|
103 nc = (m).columns (); \ |
|
104 } \ |
|
105 else \ |
|
106 abort (); |
|
107 |
|
108 // Assign a real or complex value to a tree_constant. |
|
109 // |
|
110 // This macro only makes sense inside a friend or member function of |
|
111 // the tree_constant_rep class. |
|
112 |
|
113 #define REP_ELEM_ASSIGN(i,j,rval,cval,real_type) \ |
|
114 do \ |
|
115 { \ |
|
116 if (type_tag == tree_constant_rep::matrix_constant) \ |
|
117 { \ |
|
118 if (real_type) \ |
|
119 matrix->elem ((i), (j)) = (rval); \ |
|
120 else \ |
|
121 abort (); \ |
|
122 } \ |
|
123 else \ |
|
124 { \ |
|
125 if (real_type) \ |
|
126 complex_matrix->elem ((i), (j)) = (rval); \ |
|
127 else \ |
|
128 complex_matrix->elem ((i), (j)) = (cval); \ |
|
129 } \ |
|
130 } \ |
|
131 while (0) |
|
132 |
|
133 // Given a real and complex matrix and row and column dimensions, |
|
134 // declare both and size one of them. Only one of the matrices should |
|
135 // be used after this macro has been used. |
|
136 |
|
137 // This macro only makes sense inside a friend or member function of |
|
138 // the tree_constant_rep class. |
|
139 |
|
140 #define CRMATRIX(m,cm,nr,nc) \ |
|
141 Matrix m; \ |
|
142 ComplexMatrix cm; \ |
|
143 if (type_tag == tree_constant_rep::matrix_constant) \ |
|
144 (m).resize ((nr), (nc)); \ |
|
145 else if (type_tag == complex_matrix_constant) \ |
|
146 (cm).resize ((nr), (nc)); \ |
|
147 else \ |
|
148 abort (); \ |
|
149 |
|
150 // Assign a real or complex matrix to a tree constant. |
|
151 |
|
152 // This macro only makes sense inside a friend or member function of |
|
153 // the tree_constant_rep class. |
|
154 |
|
155 #define ASSIGN_CRMATRIX_TO(tc,m,cm) \ |
|
156 do \ |
|
157 { \ |
|
158 if (type_tag == matrix_constant) \ |
|
159 tc = tree_constant (m); \ |
|
160 else \ |
|
161 tc = tree_constant (cm); \ |
|
162 } \ |
|
163 while (0) |
|
164 |
|
165 // Assign an element of this tree_constant_rep's real or complex |
|
166 // matrix to another real or complex matrix. |
|
167 |
|
168 // This macro only makes sense inside a friend or member function of |
|
169 // the tree_constant_rep class. |
|
170 |
|
171 #define CRMATRIX_ASSIGN_REP_ELEM(m,cm,i1,j1,i2,j2) \ |
|
172 do \ |
|
173 { \ |
|
174 if (type_tag == matrix_constant) \ |
|
175 (m).elem ((i1), (j1)) = matrix->elem ((i2), (j2)); \ |
|
176 else \ |
|
177 (cm).elem ((i1), (j1)) = complex_matrix->elem ((i2), (j2)); \ |
|
178 } \ |
|
179 while (0) |
|
180 |
|
181 // Assign a value to an element of a real or complex matrix. Assumes |
|
182 // that the lhs and rhs are either both real or both complex types. |
|
183 |
|
184 #define CRMATRIX_ASSIGN_ELEM(m,cm,i,j,rval,cval,real_type) \ |
|
185 do \ |
|
186 { \ |
|
187 if (real_type) \ |
|
188 (m).elem ((i), (j)) = (rval); \ |
|
189 else \ |
|
190 (cm).elem ((i), (j)) = (cval); \ |
|
191 } \ |
|
192 while (0) |
|
193 |
|
194 |
|
195 // A couple of handy helper functions. |
|
196 |
|
197 static int |
|
198 any_element_less_than (const Matrix& a, double val) |
|
199 { |
|
200 int nr = a.rows (); |
|
201 int nc = a.columns (); |
|
202 for (int j = 0; j < nc; j++) |
|
203 for (int i = 0; i < nr; i++) |
|
204 if (a.elem (i, j) < val) |
|
205 return 1; |
|
206 return 0; |
|
207 } |
|
208 |
|
209 static int |
|
210 any_element_greater_than (const Matrix& a, double val) |
|
211 { |
|
212 int nr = a.rows (); |
|
213 int nc = a.columns (); |
|
214 for (int j = 0; j < nc; j++) |
|
215 for (int i = 0; i < nr; i++) |
|
216 if (a.elem (i, j) > val) |
|
217 return 1; |
|
218 return 0; |
|
219 } |
|
220 |
|
221 static int |
|
222 any_element_is_complex (const ComplexMatrix& a) |
|
223 { |
|
224 int nr = a.rows (); |
|
225 int nc = a.columns (); |
|
226 for (int j = 0; j < nc; j++) |
|
227 for (int i = 0; i < nr; i++) |
|
228 if (imag (a.elem (i, j)) != 0.0) |
|
229 return 1; |
|
230 return 0; |
|
231 } |
|
232 |
|
233 // Now, the classes. |
|
234 |
|
235 /* |
|
236 * The real representation of constants. |
|
237 */ |
|
238 tree_constant_rep::tree_constant_rep (void) |
|
239 { |
|
240 type_tag = unknown_constant; |
|
241 } |
|
242 |
|
243 tree_constant_rep::tree_constant_rep (double d) |
|
244 { |
|
245 scalar = d; |
|
246 type_tag = scalar_constant; |
|
247 } |
|
248 |
|
249 tree_constant_rep::tree_constant_rep (const Matrix& m) |
|
250 { |
|
251 if (m.rows () == 1 && m.columns () == 1) |
|
252 { |
|
253 scalar = m.elem (0, 0); |
|
254 type_tag = scalar_constant; |
|
255 } |
|
256 else |
|
257 { |
|
258 matrix = new Matrix (m); |
|
259 type_tag = matrix_constant; |
|
260 } |
|
261 } |
|
262 |
|
263 tree_constant_rep::tree_constant_rep (const DiagMatrix& d) |
|
264 { |
|
265 if (d.rows () == 1 && d.columns () == 1) |
|
266 { |
|
267 scalar = d.elem (0, 0); |
|
268 type_tag = scalar_constant; |
|
269 } |
|
270 else |
|
271 { |
|
272 matrix = new Matrix (d); |
|
273 type_tag = matrix_constant; |
|
274 } |
|
275 } |
|
276 |
|
277 tree_constant_rep::tree_constant_rep (const RowVector& v, int |
|
278 prefer_column_vector) |
|
279 { |
|
280 int len = v.capacity (); |
|
281 if (len == 1) |
|
282 { |
|
283 scalar = v.elem (0); |
|
284 type_tag = scalar_constant; |
|
285 } |
|
286 else |
|
287 { |
|
288 int pcv = (prefer_column_vector < 0) |
|
289 ? user_pref.prefer_column_vectors |
|
290 : prefer_column_vector; |
|
291 |
|
292 if (pcv) |
|
293 { |
|
294 Matrix m (len, 1); |
|
295 for (int i = 0; i < len; i++) |
|
296 m.elem (i, 0) = v.elem (i); |
|
297 matrix = new Matrix (m); |
|
298 type_tag = matrix_constant; |
|
299 } |
|
300 else |
|
301 { |
|
302 Matrix m (1, len); |
|
303 for (int i = 0; i < len; i++) |
|
304 m.elem (0, i) = v.elem (i); |
|
305 matrix = new Matrix (m); |
|
306 type_tag = matrix_constant; |
|
307 } |
|
308 } |
|
309 } |
|
310 |
|
311 tree_constant_rep::tree_constant_rep (const ColumnVector& v, |
|
312 int prefer_column_vector) |
|
313 { |
|
314 int len = v.capacity (); |
|
315 if (len == 1) |
|
316 { |
|
317 scalar = v.elem (0); |
|
318 type_tag = scalar_constant; |
|
319 } |
|
320 else |
|
321 { |
|
322 int pcv = (prefer_column_vector < 0) |
|
323 ? user_pref.prefer_column_vectors |
|
324 : prefer_column_vector; |
|
325 |
|
326 if (pcv) |
|
327 { |
|
328 Matrix m (len, 1); |
|
329 for (int i = 0; i < len; i++) |
|
330 m.elem (i, 0) = v.elem (i); |
|
331 matrix = new Matrix (m); |
|
332 type_tag = matrix_constant; |
|
333 } |
|
334 else |
|
335 { |
|
336 Matrix m (1, len); |
|
337 for (int i = 0; i < len; i++) |
|
338 m.elem (0, i) = v.elem (i); |
|
339 matrix = new Matrix (m); |
|
340 type_tag = matrix_constant; |
|
341 } |
|
342 } |
|
343 } |
|
344 |
|
345 tree_constant_rep::tree_constant_rep (const Complex& c) |
|
346 { |
|
347 complex_scalar = new Complex (c); |
|
348 type_tag = complex_scalar_constant; |
|
349 } |
|
350 |
|
351 tree_constant_rep::tree_constant_rep (const ComplexMatrix& m) |
|
352 { |
|
353 if (m.rows () == 1 && m.columns () == 1) |
|
354 { |
|
355 complex_scalar = new Complex (m.elem (0, 0)); |
|
356 type_tag = complex_scalar_constant; |
|
357 } |
|
358 else |
|
359 { |
|
360 complex_matrix = new ComplexMatrix (m); |
|
361 type_tag = complex_matrix_constant; |
|
362 } |
|
363 } |
|
364 |
|
365 tree_constant_rep::tree_constant_rep (const ComplexDiagMatrix& d) |
|
366 { |
|
367 if (d.rows () == 1 && d.columns () == 1) |
|
368 { |
|
369 complex_scalar = new Complex (d.elem (0, 0)); |
|
370 type_tag = complex_scalar_constant; |
|
371 } |
|
372 else |
|
373 { |
|
374 complex_matrix = new ComplexMatrix (d); |
|
375 type_tag = complex_matrix_constant; |
|
376 } |
|
377 } |
|
378 |
|
379 tree_constant_rep::tree_constant_rep (const ComplexRowVector& v, |
|
380 int prefer_column_vector) |
|
381 { |
|
382 int len = v.capacity (); |
|
383 if (len == 1) |
|
384 { |
|
385 complex_scalar = new Complex (v.elem (0)); |
|
386 type_tag = complex_scalar_constant; |
|
387 } |
|
388 else |
|
389 { |
|
390 int pcv = (prefer_column_vector < 0) |
|
391 ? user_pref.prefer_column_vectors |
|
392 : prefer_column_vector; |
|
393 |
|
394 if (pcv) |
|
395 { |
|
396 ComplexMatrix m (len, 1); |
|
397 for (int i = 0; i < len; i++) |
|
398 m.elem (i, 0) = v.elem (i); |
|
399 complex_matrix = new ComplexMatrix (m); |
|
400 type_tag = complex_matrix_constant; |
|
401 } |
|
402 else |
|
403 { |
|
404 ComplexMatrix m (1, len); |
|
405 for (int i = 0; i < len; i++) |
|
406 m.elem (0, i) = v.elem (i); |
|
407 complex_matrix = new ComplexMatrix (m); |
|
408 type_tag = complex_matrix_constant; |
|
409 } |
|
410 } |
|
411 } |
|
412 |
|
413 tree_constant_rep::tree_constant_rep (const ComplexColumnVector& v, |
|
414 int prefer_column_vector) |
|
415 { |
|
416 int len = v.capacity (); |
|
417 if (len == 1) |
|
418 { |
|
419 complex_scalar = new Complex (v.elem (0)); |
|
420 type_tag = complex_scalar_constant; |
|
421 } |
|
422 else |
|
423 { |
|
424 int pcv = (prefer_column_vector < 0) |
|
425 ? user_pref.prefer_column_vectors |
|
426 : prefer_column_vector; |
|
427 |
|
428 if (pcv) |
|
429 { |
|
430 ComplexMatrix m (len, 1); |
|
431 for (int i = 0; i < len; i++) |
|
432 m.elem (i, 0) = v.elem (i); |
|
433 complex_matrix = new ComplexMatrix (m); |
|
434 type_tag = complex_matrix_constant; |
|
435 } |
|
436 else |
|
437 { |
|
438 ComplexMatrix m (1, len); |
|
439 for (int i = 0; i < len; i++) |
|
440 m.elem (0, i) = v.elem (i); |
|
441 complex_matrix = new ComplexMatrix (m); |
|
442 type_tag = complex_matrix_constant; |
|
443 } |
|
444 } |
|
445 } |
|
446 |
|
447 tree_constant_rep::tree_constant_rep (const char *s) |
|
448 { |
|
449 string = strsave (s); |
|
450 type_tag = string_constant; |
|
451 } |
|
452 |
|
453 tree_constant_rep::tree_constant_rep (double b, double l, double i) |
|
454 { |
|
455 range = new Range (b, l, i); |
|
456 int nel = range->nelem (); |
|
457 if (nel < 0) |
|
458 { |
|
459 delete range; |
|
460 type_tag = unknown_constant; |
|
461 if (nel == -1) |
|
462 ::error ("number of elements in range exceeds INT_MAX"); |
|
463 else |
|
464 ::error ("invalid range"); |
|
465 } |
|
466 else if (nel > 1) |
|
467 type_tag = range_constant; |
|
468 else |
|
469 { |
|
470 delete range; |
|
471 if (nel == 1) |
|
472 { |
|
473 scalar = b; |
|
474 type_tag = scalar_constant; |
|
475 } |
|
476 else if (nel == 0) |
|
477 { |
|
478 matrix = new Matrix (); |
|
479 type_tag = matrix_constant; |
|
480 } |
|
481 else |
|
482 panic_impossible (); |
|
483 } |
|
484 } |
|
485 |
|
486 tree_constant_rep::tree_constant_rep (const Range& r) |
|
487 { |
|
488 if (r.nelem () > 1) |
|
489 { |
|
490 range = new Range (r); |
|
491 type_tag = range_constant; |
|
492 } |
|
493 else if (r.nelem () == 1) |
|
494 { |
|
495 scalar = r.base (); |
|
496 type_tag = scalar_constant; |
|
497 } |
|
498 else if (r.nelem () == 0) |
|
499 { |
|
500 matrix = new Matrix (); |
|
501 type_tag = matrix_constant; |
|
502 } |
|
503 else |
|
504 panic_impossible (); |
|
505 } |
|
506 |
|
507 tree_constant_rep::tree_constant_rep (tree_constant_rep::constant_type t) |
|
508 { |
|
509 assert (t == magic_colon); |
|
510 |
|
511 type_tag = magic_colon; |
|
512 } |
|
513 |
|
514 tree_constant_rep::tree_constant_rep (const tree_constant_rep& t) |
|
515 { |
|
516 type_tag = t.type_tag; |
|
517 |
|
518 switch (t.type_tag) |
|
519 { |
|
520 case unknown_constant: |
|
521 break; |
|
522 case scalar_constant: |
|
523 scalar = t.scalar; |
|
524 break; |
|
525 case matrix_constant: |
|
526 matrix = new Matrix (*(t.matrix)); |
|
527 break; |
|
528 case string_constant: |
|
529 string = strsave (t.string); |
|
530 break; |
|
531 case complex_matrix_constant: |
|
532 complex_matrix = new ComplexMatrix (*(t.complex_matrix)); |
|
533 break; |
|
534 case complex_scalar_constant: |
|
535 complex_scalar = new Complex (*(t.complex_scalar)); |
|
536 break; |
|
537 case range_constant: |
|
538 range = new Range (*(t.range)); |
|
539 break; |
|
540 case magic_colon: |
|
541 break; |
|
542 default: |
|
543 panic_impossible (); |
|
544 break; |
|
545 } |
|
546 } |
|
547 |
|
548 tree_constant_rep::~tree_constant_rep (void) |
|
549 { |
|
550 switch (type_tag) |
|
551 { |
|
552 case unknown_constant: |
|
553 break; |
|
554 case scalar_constant: |
|
555 break; |
|
556 case matrix_constant: |
|
557 delete matrix; |
|
558 break; |
|
559 case complex_scalar_constant: |
|
560 delete complex_scalar; |
|
561 break; |
|
562 case complex_matrix_constant: |
|
563 delete complex_matrix; |
|
564 break; |
|
565 case string_constant: |
|
566 delete [] string; |
|
567 break; |
|
568 case range_constant: |
|
569 delete range; |
|
570 break; |
|
571 case magic_colon: |
|
572 break; |
|
573 default: |
|
574 panic_impossible (); |
|
575 break; |
|
576 } |
|
577 } |
|
578 |
|
579 #if defined (MDEBUG) |
|
580 void * |
|
581 tree_constant_rep::operator new (size_t size) |
|
582 { |
|
583 tree_constant_rep *p = ::new tree_constant_rep; |
|
584 cerr << "tree_constant_rep::new(): " << p << "\n"; |
|
585 return p; |
|
586 } |
|
587 |
|
588 void |
|
589 tree_constant_rep::operator delete (void *p, size_t size) |
|
590 { |
|
591 cerr << "tree_constant_rep::delete(): " << p << "\n"; |
|
592 ::delete p; |
|
593 } |
|
594 #endif |
|
595 |
|
596 void |
|
597 tree_constant_rep::resize (int i, int j) |
|
598 { |
|
599 switch (type_tag) |
|
600 { |
|
601 case matrix_constant: |
|
602 matrix->resize (i, j); |
|
603 break; |
|
604 case complex_matrix_constant: |
|
605 complex_matrix->resize (i, j); |
|
606 break; |
|
607 default: |
|
608 panic_impossible (); |
|
609 break; |
|
610 } |
|
611 } |
|
612 |
|
613 void |
|
614 tree_constant_rep::resize (int i, int j, double val) |
|
615 { |
|
616 switch (type_tag) |
|
617 { |
|
618 case matrix_constant: |
|
619 matrix->resize (i, j, val); |
|
620 break; |
|
621 case complex_matrix_constant: |
|
622 complex_matrix->resize (i, j, val); |
|
623 break; |
|
624 default: |
|
625 panic_impossible (); |
|
626 break; |
|
627 } |
|
628 } |
|
629 |
|
630 void |
|
631 tree_constant_rep::maybe_resize (int i, int j) |
|
632 { |
|
633 int nr = rows (); |
|
634 int nc = columns (); |
|
635 |
|
636 i++; |
|
637 j++; |
|
638 |
|
639 assert (i > 0 && j > 0); |
|
640 |
|
641 if (i > nr || j > nc) |
|
642 { |
|
643 if (user_pref.resize_on_range_error) |
|
644 resize (MAX (i, nr), MAX (j, nc), 0.0); |
|
645 else |
|
646 { |
|
647 if (i > nr) |
|
648 ::error ("row index = %d exceeds max row dimension = %d", i, nr); |
|
649 |
|
650 if (j > nc) |
|
651 ::error ("column index = %d exceeds max column dimension = %d", |
|
652 j, nc); |
|
653 } |
|
654 } |
|
655 } |
|
656 |
|
657 void |
|
658 tree_constant_rep::maybe_resize (int i, force_orient f_orient = no_orient) |
|
659 { |
|
660 int nr = rows (); |
|
661 int nc = columns (); |
|
662 |
|
663 i++; |
|
664 |
|
665 assert (i >= 0 && (nr <= 1 || nc <= 1)); |
|
666 |
|
667 // This function never reduces the size of a vector, and all vectors |
|
668 // have dimensions of at least 0x0. If i is 0, it is either because |
|
669 // a vector has been indexed with a vector of all zeros (in which case |
|
670 // the index vector is empty and nothing will happen) or a vector has |
|
671 // been indexed with 0 (an error which will be caught elsewhere). |
|
672 if (i == 0) |
|
673 return; |
|
674 |
|
675 if (nr <= 1 && nc <= 1 && i >= 1) |
|
676 { |
|
677 if (user_pref.resize_on_range_error) |
|
678 { |
|
679 if (f_orient == row_orient) |
|
680 resize (1, i, 0.0); |
|
681 else if (f_orient == column_orient) |
|
682 resize (i, 1, 0.0); |
|
683 else if (user_pref.prefer_column_vectors) |
|
684 resize (i, 1, 0.0); |
|
685 else |
|
686 resize (1, i, 0.0); |
|
687 } |
|
688 else |
|
689 ::error ("matrix index = %d exceeds max dimension = %d", i, nc); |
|
690 } |
|
691 else if (nr == 1 && i > nc) |
|
692 { |
|
693 if (user_pref.resize_on_range_error) |
|
694 resize (1, i, 0.0); |
|
695 else |
|
696 ::error ("matrix index = %d exceeds max dimension = %d", i, nc); |
|
697 } |
|
698 else if (nc == 1 && i > nr) |
|
699 { |
|
700 if (user_pref.resize_on_range_error) |
|
701 resize (i, 1, 0.0); |
|
702 else |
|
703 ::error ("matrix index = %d exceeds max dimension = ", i, nc); |
|
704 } |
|
705 } |
|
706 |
|
707 double |
|
708 tree_constant_rep::to_scalar (void) const |
|
709 { |
|
710 tree_constant tmp = make_numeric (); |
|
711 |
|
712 double retval = 0.0; |
|
713 |
|
714 switch (tmp.const_type ()) |
|
715 { |
|
716 case tree_constant_rep::scalar_constant: |
|
717 case tree_constant_rep::complex_scalar_constant: |
|
718 retval = tmp.double_value (); |
|
719 break; |
|
720 case tree_constant_rep::matrix_constant: |
|
721 if (user_pref.do_fortran_indexing) |
|
722 { |
|
723 Matrix m = tmp.matrix_value (); |
|
724 retval = m (0, 0); |
|
725 } |
|
726 break; |
|
727 case tree_constant_rep::complex_matrix_constant: |
|
728 if (user_pref.do_fortran_indexing) |
|
729 { |
|
730 int flag = user_pref.ok_to_lose_imaginary_part; |
|
731 if (flag == -1) |
|
732 warning ("implicit conversion of complex value to real value"); |
|
733 |
|
734 if (flag != 0) |
|
735 { |
|
736 ComplexMatrix m = tmp.complex_matrix_value (); |
|
737 return ::real (m (0, 0)); |
|
738 } |
|
739 else |
|
740 jump_to_top_level (); |
|
741 } |
|
742 else |
|
743 { |
|
744 ::error ("complex matrix used in invalid context"); |
|
745 jump_to_top_level (); |
|
746 } |
|
747 break; |
|
748 default: |
|
749 break; |
|
750 } |
|
751 return retval; |
|
752 } |
|
753 |
|
754 ColumnVector |
|
755 tree_constant_rep::to_vector (void) const |
|
756 { |
|
757 tree_constant tmp = make_numeric (); |
|
758 |
|
759 ColumnVector retval; |
|
760 |
|
761 switch (tmp.const_type ()) |
|
762 { |
|
763 case tree_constant_rep::scalar_constant: |
|
764 case tree_constant_rep::complex_scalar_constant: |
|
765 retval.resize (1); |
|
766 retval.elem (0) = tmp.double_value (); |
|
767 break; |
|
768 case tree_constant_rep::complex_matrix_constant: |
|
769 case tree_constant_rep::matrix_constant: |
|
770 { |
|
771 Matrix m = tmp.matrix_value (); |
|
772 int nr = m.rows (); |
|
773 int nc = m.columns (); |
|
774 if (nr == 1) |
|
775 { |
|
776 retval.resize (nc); |
|
777 for (int i = 0; i < nc; i++) |
|
778 retval.elem (i) = m (0, i); |
|
779 } |
|
780 else if (nc == 1) |
|
781 { |
|
782 retval.resize (nr); |
|
783 for (int i = 0; i < nr; i++) |
|
784 retval.elem (i) = m.elem (i, 0); |
|
785 } |
|
786 } |
|
787 break; |
|
788 default: |
|
789 panic_impossible (); |
|
790 break; |
|
791 } |
|
792 return retval; |
|
793 } |
|
794 |
|
795 Matrix |
|
796 tree_constant_rep::to_matrix (void) const |
|
797 { |
|
798 tree_constant tmp = make_numeric (); |
|
799 |
|
800 Matrix retval; |
|
801 |
|
802 switch (tmp.const_type ()) |
|
803 { |
|
804 case tree_constant_rep::scalar_constant: |
|
805 retval.resize (1, 1); |
|
806 retval.elem (0, 0) = tmp.double_value (); |
|
807 break; |
|
808 case tree_constant_rep::matrix_constant: |
|
809 retval = tmp.matrix_value (); |
|
810 break; |
|
811 default: |
|
812 break; |
|
813 } |
|
814 return retval; |
|
815 } |
|
816 |
|
817 tree_constant_rep::constant_type |
|
818 tree_constant_rep::force_numeric (int force_str_conv = 0) |
|
819 { |
|
820 switch (type_tag) |
|
821 { |
|
822 case scalar_constant: |
|
823 case matrix_constant: |
|
824 case complex_scalar_constant: |
|
825 case complex_matrix_constant: |
|
826 break; |
|
827 case string_constant: |
|
828 { |
|
829 if (! force_str_conv && ! user_pref.implicit_str_to_num_ok) |
|
830 { |
|
831 ::error ("failed to convert `%s' to a numeric type --", string); |
|
832 ::error ("default conversion turned off"); |
|
833 // Abort! |
|
834 jump_to_top_level (); |
|
835 } |
|
836 |
|
837 int len = strlen (string); |
|
838 if (len > 1) |
|
839 { |
|
840 type_tag = matrix_constant; |
|
841 Matrix *tm = new Matrix (1, len); |
|
842 for (int i = 0; i < len; i++) |
|
843 tm->elem (0, i) = toascii ((int) string[i]); |
|
844 matrix = tm; |
|
845 } |
|
846 else if (len == 1) |
|
847 { |
|
848 type_tag = scalar_constant; |
|
849 scalar = toascii ((int) string[0]); |
|
850 } |
|
851 else if (len == 0) |
|
852 { |
|
853 type_tag = matrix_constant; |
|
854 matrix = new Matrix (0, 0); |
|
855 } |
|
856 else |
|
857 panic_impossible (); |
|
858 } |
|
859 break; |
|
860 case range_constant: |
|
861 { |
|
862 int len = range->nelem (); |
|
863 if (len > 1) |
|
864 { |
|
865 type_tag = matrix_constant; |
|
866 Matrix *tm = new Matrix (1, len); |
|
867 double b = range->base (); |
|
868 double increment = range->inc (); |
|
869 for (int i = 0; i < len; i++) |
|
870 tm->elem (0, i) = b + i * increment; |
|
871 matrix = tm; |
|
872 } |
|
873 else if (len == 1) |
|
874 { |
|
875 type_tag = scalar_constant; |
|
876 scalar = range->base (); |
|
877 } |
|
878 } |
|
879 break; |
|
880 case magic_colon: |
|
881 default: |
|
882 panic_impossible (); |
|
883 break; |
|
884 } |
|
885 return type_tag; |
|
886 } |
|
887 |
|
888 tree_constant |
|
889 tree_constant_rep::make_numeric (int force_str_conv = 0) const |
|
890 { |
|
891 tree_constant retval; |
|
892 switch (type_tag) |
|
893 { |
|
894 case scalar_constant: |
|
895 retval = tree_constant (scalar); |
|
896 break; |
|
897 case matrix_constant: |
|
898 retval = tree_constant (*matrix); |
|
899 break; |
|
900 case complex_scalar_constant: |
|
901 retval = tree_constant (*complex_scalar); |
|
902 break; |
|
903 case complex_matrix_constant: |
|
904 retval = tree_constant (*complex_matrix); |
|
905 break; |
|
906 case string_constant: |
|
907 retval = tree_constant (string); |
|
908 retval.force_numeric (force_str_conv); |
|
909 break; |
|
910 case range_constant: |
|
911 retval = tree_constant (*range); |
|
912 retval.force_numeric (force_str_conv); |
|
913 break; |
|
914 case magic_colon: |
|
915 default: |
|
916 panic_impossible (); |
|
917 break; |
|
918 } |
|
919 return retval; |
|
920 } |
|
921 |
|
922 tree_constant |
|
923 do_binary_op (tree_constant& a, tree_constant& b, tree::expression_type t) |
|
924 { |
|
925 tree_constant ans; |
|
926 |
|
927 int first_empty = (a.rows () == 0 || a.columns () == 0); |
|
928 int second_empty = (b.rows () == 0 || b.columns () == 0); |
|
929 |
|
930 if (first_empty || second_empty) |
|
931 { |
|
932 int flag = user_pref.propagate_empty_matrices; |
|
933 if (flag < 0) |
|
934 warning ("binary operation on empty matrix"); |
|
935 else if (flag == 0) |
|
936 { |
|
937 ::error ("invalid binary operation on empty matrix"); |
|
938 return ans; |
|
939 } |
|
940 } |
|
941 |
|
942 tree_constant tmp_a = a.make_numeric (); |
|
943 tree_constant tmp_b = b.make_numeric (); |
|
944 |
|
945 tree_constant_rep::constant_type a_type = tmp_a.const_type (); |
|
946 tree_constant_rep::constant_type b_type = tmp_b.const_type (); |
|
947 |
|
948 double d1, d2; |
|
949 Matrix m1, m2; |
|
950 Complex c1, c2; |
|
951 ComplexMatrix cm1, cm2; |
|
952 |
|
953 switch (a_type) |
|
954 { |
|
955 case tree_constant_rep::scalar_constant: |
|
956 d1 = tmp_a.double_value (); |
|
957 switch (b_type) |
|
958 { |
|
959 case tree_constant_rep::scalar_constant: |
|
960 d2 = tmp_b.double_value (); |
|
961 ans = do_binary_op (d1, d2, t); |
|
962 break; |
|
963 case tree_constant_rep::matrix_constant: |
|
964 m2 = tmp_b.matrix_value (); |
|
965 ans = do_binary_op (d1, m2, t); |
|
966 break; |
|
967 case tree_constant_rep::complex_scalar_constant: |
|
968 c2 = tmp_b.complex_value (); |
|
969 ans = do_binary_op (d1, c2, t); |
|
970 break; |
|
971 case tree_constant_rep::complex_matrix_constant: |
|
972 cm2 = tmp_b.complex_matrix_value (); |
|
973 ans = do_binary_op (d1, cm2, t); |
|
974 break; |
|
975 case tree_constant_rep::magic_colon: |
|
976 default: |
|
977 panic_impossible (); |
|
978 break; |
|
979 } |
|
980 break; |
|
981 case tree_constant_rep::matrix_constant: |
|
982 m1 = tmp_a.matrix_value (); |
|
983 switch (b_type) |
|
984 { |
|
985 case tree_constant_rep::scalar_constant: |
|
986 d2 = tmp_b.double_value (); |
|
987 ans = do_binary_op (m1, d2, t); |
|
988 break; |
|
989 case tree_constant_rep::matrix_constant: |
|
990 m2 = tmp_b.matrix_value (); |
|
991 ans = do_binary_op (m1, m2, t); |
|
992 break; |
|
993 case tree_constant_rep::complex_scalar_constant: |
|
994 c2 = tmp_b.complex_value (); |
|
995 ans = do_binary_op (m1, c2, t); |
|
996 break; |
|
997 case tree_constant_rep::complex_matrix_constant: |
|
998 cm2 = tmp_b.complex_matrix_value (); |
|
999 ans = do_binary_op (m1, cm2, t); |
|
1000 break; |
|
1001 case tree_constant_rep::magic_colon: |
|
1002 default: |
|
1003 panic_impossible (); |
|
1004 break; |
|
1005 } |
|
1006 break; |
|
1007 case tree_constant_rep::complex_scalar_constant: |
|
1008 c1 = tmp_a.complex_value (); |
|
1009 switch (b_type) |
|
1010 { |
|
1011 case tree_constant_rep::scalar_constant: |
|
1012 d2 = tmp_b.double_value (); |
|
1013 ans = do_binary_op (c1, d2, t); |
|
1014 break; |
|
1015 case tree_constant_rep::matrix_constant: |
|
1016 m2 = tmp_b.matrix_value (); |
|
1017 ans = do_binary_op (c1, m2, t); |
|
1018 break; |
|
1019 case tree_constant_rep::complex_scalar_constant: |
|
1020 c2 = tmp_b.complex_value (); |
|
1021 ans = do_binary_op (c1, c2, t); |
|
1022 break; |
|
1023 case tree_constant_rep::complex_matrix_constant: |
|
1024 cm2 = tmp_b.complex_matrix_value (); |
|
1025 ans = do_binary_op (c1, cm2, t); |
|
1026 break; |
|
1027 case tree_constant_rep::magic_colon: |
|
1028 default: |
|
1029 panic_impossible (); |
|
1030 break; |
|
1031 } |
|
1032 break; |
|
1033 case tree_constant_rep::complex_matrix_constant: |
|
1034 cm1 = tmp_a.complex_matrix_value (); |
|
1035 switch (b_type) |
|
1036 { |
|
1037 case tree_constant_rep::scalar_constant: |
|
1038 d2 = tmp_b.double_value (); |
|
1039 ans = do_binary_op (cm1, d2, t); |
|
1040 break; |
|
1041 case tree_constant_rep::matrix_constant: |
|
1042 m2 = tmp_b.matrix_value (); |
|
1043 ans = do_binary_op (cm1, m2, t); |
|
1044 break; |
|
1045 case tree_constant_rep::complex_scalar_constant: |
|
1046 c2 = tmp_b.complex_value (); |
|
1047 ans = do_binary_op (cm1, c2, t); |
|
1048 break; |
|
1049 case tree_constant_rep::complex_matrix_constant: |
|
1050 cm2 = tmp_b.complex_matrix_value (); |
|
1051 ans = do_binary_op (cm1, cm2, t); |
|
1052 break; |
|
1053 case tree_constant_rep::magic_colon: |
|
1054 default: |
|
1055 panic_impossible (); |
|
1056 break; |
|
1057 } |
|
1058 break; |
|
1059 case tree_constant_rep::magic_colon: |
|
1060 default: |
|
1061 panic_impossible (); |
|
1062 break; |
|
1063 } |
|
1064 return ans; |
|
1065 } |
|
1066 |
|
1067 tree_constant |
|
1068 do_unary_op (tree_constant& a, tree::expression_type t) |
|
1069 { |
|
1070 tree_constant ans; |
|
1071 |
|
1072 if (a.rows () == 0 || a.columns () == 0) |
|
1073 { |
|
1074 int flag = user_pref.propagate_empty_matrices; |
|
1075 if (flag < 0) |
|
1076 warning ("unary operation on empty matrix"); |
|
1077 else if (flag == 0) |
|
1078 { |
|
1079 ::error ("invalid unary operation on empty matrix"); |
|
1080 return ans; |
|
1081 } |
|
1082 } |
|
1083 |
|
1084 tree_constant tmp_a = a.make_numeric (); |
|
1085 |
|
1086 switch (tmp_a.const_type ()) |
|
1087 { |
|
1088 case tree_constant_rep::scalar_constant: |
|
1089 ans = do_unary_op (tmp_a.double_value (), t); |
|
1090 break; |
|
1091 case tree_constant_rep::matrix_constant: |
|
1092 { |
|
1093 Matrix m = tmp_a.matrix_value (); |
|
1094 ans = do_unary_op (m, t); |
|
1095 } |
|
1096 break; |
|
1097 case tree_constant_rep::complex_scalar_constant: |
|
1098 ans = do_unary_op (tmp_a.complex_value (), t); |
|
1099 break; |
|
1100 case tree_constant_rep::complex_matrix_constant: |
|
1101 { |
|
1102 ComplexMatrix m = tmp_a.complex_matrix_value (); |
|
1103 ans = do_unary_op (m, t); |
|
1104 } |
|
1105 break; |
|
1106 case tree_constant_rep::magic_colon: |
|
1107 default: |
|
1108 panic_impossible (); |
|
1109 break; |
|
1110 } |
|
1111 return ans; |
|
1112 } |
|
1113 |
|
1114 void |
|
1115 tree_constant_rep::bump_value (tree::expression_type etype) |
|
1116 { |
|
1117 switch (etype) |
|
1118 { |
|
1119 case tree::increment: |
|
1120 switch (type_tag) |
|
1121 { |
|
1122 case scalar_constant: |
|
1123 scalar++; |
|
1124 break; |
|
1125 case matrix_constant: |
|
1126 *matrix = *matrix + 1.0; |
|
1127 break; |
|
1128 case complex_scalar_constant: |
|
1129 *complex_scalar = *complex_scalar + 1.0; |
|
1130 break; |
|
1131 case complex_matrix_constant: |
|
1132 *complex_matrix = *complex_matrix + 1.0; |
|
1133 break; |
|
1134 case string_constant: |
|
1135 ::error ("string++ and ++string not implemented yet, ok?"); |
|
1136 break; |
|
1137 case range_constant: |
|
1138 range->set_base (range->base () + 1.0); |
|
1139 range->set_limit (range->limit () + 1.0); |
|
1140 break; |
|
1141 case magic_colon: |
|
1142 default: |
|
1143 panic_impossible (); |
|
1144 break; |
|
1145 } |
|
1146 break; |
|
1147 case tree::decrement: |
|
1148 switch (type_tag) |
|
1149 { |
|
1150 case scalar_constant: |
|
1151 scalar--; |
|
1152 break; |
|
1153 case matrix_constant: |
|
1154 *matrix = *matrix - 1.0; |
|
1155 break; |
|
1156 case string_constant: |
|
1157 ::error ("string-- and -- string not implemented yet, ok?"); |
|
1158 break; |
|
1159 case range_constant: |
|
1160 range->set_base (range->base () - 1.0); |
|
1161 range->set_limit (range->limit () - 1.0); |
|
1162 break; |
|
1163 case magic_colon: |
|
1164 default: |
|
1165 panic_impossible (); |
|
1166 break; |
|
1167 } |
|
1168 break; |
|
1169 default: |
|
1170 panic_impossible (); |
|
1171 break; |
|
1172 } |
|
1173 } |
|
1174 |
|
1175 void |
|
1176 tree_constant_rep::maybe_mutate (void) |
|
1177 { |
|
1178 if (error_state) |
|
1179 return; |
|
1180 |
|
1181 switch (type_tag) |
|
1182 { |
|
1183 case complex_scalar_constant: |
|
1184 if (::imag (*complex_scalar) == 0.0) |
|
1185 { |
|
1186 double d = ::real (*complex_scalar); |
|
1187 delete complex_scalar; |
|
1188 scalar = d; |
|
1189 type_tag = scalar_constant; |
|
1190 } |
|
1191 break; |
|
1192 case complex_matrix_constant: |
|
1193 if (! any_element_is_complex (*complex_matrix)) |
|
1194 { |
|
1195 Matrix *m = new Matrix (::real (*complex_matrix)); |
|
1196 delete complex_matrix; |
|
1197 matrix = m; |
|
1198 type_tag = matrix_constant; |
|
1199 } |
|
1200 break; |
|
1201 case scalar_constant: |
|
1202 case matrix_constant: |
|
1203 case string_constant: |
|
1204 case range_constant: |
|
1205 case magic_colon: |
|
1206 break; |
|
1207 default: |
|
1208 panic_impossible (); |
|
1209 break; |
|
1210 } |
|
1211 |
|
1212 // Avoid calling rows() and columns() for things like magic_colon. |
|
1213 |
|
1214 int nr = 1; |
|
1215 int nc = 1; |
|
1216 if (type_tag == matrix_constant |
|
1217 || type_tag == complex_matrix_constant |
|
1218 || type_tag == range_constant) |
|
1219 { |
|
1220 nr = rows (); |
|
1221 nc = columns (); |
|
1222 } |
|
1223 |
|
1224 switch (type_tag) |
|
1225 { |
|
1226 case matrix_constant: |
|
1227 if (nr == 1 && nc == 1) |
|
1228 { |
|
1229 double d = matrix->elem (0, 0); |
|
1230 delete matrix; |
|
1231 scalar = d; |
|
1232 type_tag = scalar_constant; |
|
1233 } |
|
1234 break; |
|
1235 case complex_matrix_constant: |
|
1236 if (nr == 1 && nc == 1) |
|
1237 { |
|
1238 Complex c = complex_matrix->elem (0, 0); |
|
1239 delete complex_matrix; |
|
1240 complex_scalar = new Complex (c); |
|
1241 type_tag = complex_scalar_constant; |
|
1242 } |
|
1243 break; |
|
1244 case range_constant: |
|
1245 if (nr == 1 && nc == 1) |
|
1246 { |
|
1247 double d = range->base (); |
|
1248 delete range; |
|
1249 scalar = d; |
|
1250 type_tag = scalar_constant; |
|
1251 } |
|
1252 break; |
|
1253 default: |
|
1254 break; |
|
1255 } |
|
1256 } |
|
1257 |
|
1258 void |
|
1259 tree_constant_rep::print (void) |
|
1260 { |
|
1261 if (error_state) |
|
1262 return; |
|
1263 |
|
1264 int nr = rows (); |
|
1265 int nc = columns (); |
|
1266 |
|
1267 if (print) |
|
1268 { |
|
1269 ostrstream output_buf; |
|
1270 switch (type_tag) |
|
1271 { |
|
1272 case scalar_constant: |
|
1273 octave_print_internal (output_buf, scalar); |
|
1274 break; |
|
1275 case matrix_constant: |
|
1276 if (nr == 0 || nc == 0) |
|
1277 { |
|
1278 output_buf << "[]"; |
|
1279 if (user_pref.print_empty_dimensions) |
|
1280 output_buf << "(" << nr << "x" << nc << ")"; |
|
1281 output_buf << "\n"; |
|
1282 } |
|
1283 else |
|
1284 octave_print_internal (output_buf, *matrix); |
|
1285 break; |
|
1286 case complex_scalar_constant: |
|
1287 octave_print_internal (output_buf, *complex_scalar); |
|
1288 break; |
|
1289 case complex_matrix_constant: |
|
1290 if (nr == 0 || nc == 0) |
|
1291 { |
|
1292 output_buf << "[]"; |
|
1293 if (user_pref.print_empty_dimensions) |
|
1294 output_buf << "(" << nr << "x" << nc << ")"; |
|
1295 output_buf << "\n"; |
|
1296 } |
|
1297 else |
|
1298 octave_print_internal (output_buf, *complex_matrix); |
|
1299 break; |
|
1300 case string_constant: |
|
1301 output_buf << string << "\n"; |
|
1302 break; |
|
1303 case range_constant: |
|
1304 octave_print_internal (output_buf, *range); |
|
1305 break; |
|
1306 case magic_colon: |
|
1307 default: |
|
1308 panic_impossible (); |
|
1309 break; |
|
1310 } |
|
1311 |
|
1312 output_buf << ends; |
|
1313 maybe_page_output (output_buf); |
|
1314 } |
|
1315 } |
|
1316 |
|
1317 tree_constant |
500
|
1318 tree_constant_rep::do_index (const Octave_object& args, int nargin) |
492
|
1319 { |
|
1320 tree_constant retval; |
|
1321 |
|
1322 if (error_state) |
|
1323 return retval; |
|
1324 |
|
1325 if (rows () == 0 || columns () == 0) |
|
1326 { |
|
1327 ::error ("attempt to index empty matrix"); |
|
1328 return retval; |
|
1329 } |
|
1330 |
|
1331 switch (type_tag) |
|
1332 { |
|
1333 case complex_scalar_constant: |
|
1334 case scalar_constant: |
|
1335 retval = do_scalar_index (args, nargin); |
|
1336 break; |
|
1337 case complex_matrix_constant: |
|
1338 case matrix_constant: |
|
1339 retval = do_matrix_index (args, nargin); |
|
1340 break; |
|
1341 case string_constant: |
|
1342 gripe_string_invalid (); |
|
1343 // retval = do_string_index (args, nargin); |
|
1344 break; |
|
1345 case magic_colon: |
|
1346 case range_constant: |
|
1347 // This isn\'t great, but it\'s easier than implementing a lot of |
|
1348 // range indexing functions. |
|
1349 force_numeric (); |
|
1350 assert (type_tag != magic_colon && type_tag != range_constant); |
|
1351 retval = do_index (args, nargin); |
|
1352 break; |
|
1353 default: |
|
1354 panic_impossible (); |
|
1355 break; |
|
1356 } |
|
1357 |
|
1358 return retval; |
|
1359 } |
|
1360 |
|
1361 int |
|
1362 tree_constant_rep::save (ostream& os, int mark_as_global, int precision) |
|
1363 { |
|
1364 switch (type_tag) |
|
1365 { |
|
1366 case scalar_constant: |
|
1367 case matrix_constant: |
|
1368 case complex_scalar_constant: |
|
1369 case complex_matrix_constant: |
|
1370 case string_constant: |
|
1371 case range_constant: |
|
1372 if (mark_as_global) |
|
1373 os << "# type: global "; |
|
1374 else |
|
1375 os << "# type: "; |
|
1376 break; |
|
1377 case magic_colon: |
|
1378 default: |
|
1379 break; |
|
1380 } |
|
1381 |
|
1382 long old_precision = os.precision (); |
|
1383 os.precision (precision); |
|
1384 |
|
1385 switch (type_tag) |
|
1386 { |
|
1387 case scalar_constant: |
|
1388 os << "scalar\n" |
|
1389 << scalar << "\n"; |
|
1390 break; |
|
1391 case matrix_constant: |
|
1392 os << "matrix\n" |
|
1393 << "# rows: " << rows () << "\n" |
|
1394 << "# columns: " << columns () << "\n" |
|
1395 << *matrix ; |
|
1396 break; |
|
1397 case complex_scalar_constant: |
|
1398 os << "complex scalar\n" |
|
1399 << *complex_scalar << "\n"; |
|
1400 break; |
|
1401 case complex_matrix_constant: |
|
1402 os << "complex matrix\n" |
|
1403 << "# rows: " << rows () << "\n" |
|
1404 << "# columns: " << columns () << "\n" |
|
1405 << *complex_matrix ; |
|
1406 break; |
|
1407 case string_constant: |
|
1408 os << "string\n" |
|
1409 << "# length: " << strlen (string) << "\n" |
|
1410 << string << "\n"; |
|
1411 break; |
|
1412 case range_constant: |
|
1413 { |
|
1414 os << "range\n" |
|
1415 << "# base, limit, increment\n" |
|
1416 << range->base () << " " |
|
1417 << range->limit () << " " |
|
1418 << range->inc () << "\n"; |
|
1419 } |
|
1420 break; |
|
1421 case magic_colon: |
|
1422 default: |
|
1423 panic_impossible (); |
|
1424 break; |
|
1425 } |
|
1426 |
|
1427 os.precision (old_precision); |
|
1428 |
|
1429 // Really want to return 1 only if write is successful. |
|
1430 return 1; |
|
1431 } |
|
1432 |
|
1433 int |
|
1434 tree_constant_rep::save_three_d (ostream& os, int parametric) |
|
1435 { |
|
1436 int nr = rows (); |
|
1437 int nc = columns (); |
|
1438 |
|
1439 switch (type_tag) |
|
1440 { |
|
1441 case matrix_constant: |
|
1442 os << "# 3D data...\n" |
|
1443 << "# type: matrix\n" |
|
1444 << "# total rows: " << nr << "\n" |
|
1445 << "# total columns: " << nc << "\n"; |
|
1446 |
|
1447 if (parametric) |
|
1448 { |
|
1449 int extras = nc % 3; |
|
1450 if (extras) |
|
1451 warning ("ignoring last %d columns", extras); |
|
1452 |
|
1453 for (int i = 0; i < nc-extras; i += 3) |
|
1454 { |
|
1455 os << matrix->extract (0, i, nr-1, i+2); |
|
1456 if (i+3 < nc-extras) |
|
1457 os << "\n"; |
|
1458 } |
|
1459 } |
|
1460 else |
|
1461 { |
|
1462 for (int i = 0; i < nc; i++) |
|
1463 { |
|
1464 os << matrix->extract (0, i, nr-1, i); |
|
1465 if (i+1 < nc) |
|
1466 os << "\n"; |
|
1467 } |
|
1468 } |
|
1469 break; |
|
1470 default: |
|
1471 ::error ("for now, I can only save real matrices in 3D format"); |
|
1472 return 0; |
|
1473 break; |
|
1474 } |
|
1475 // Really want to return 1 only if write is successful. |
|
1476 return 1; |
|
1477 } |
|
1478 |
|
1479 int |
|
1480 tree_constant_rep::load (istream& is) |
|
1481 { |
|
1482 int is_global = 0; |
|
1483 |
|
1484 type_tag = unknown_constant; |
|
1485 |
|
1486 // Look for type keyword |
|
1487 |
|
1488 char *tag = extract_keyword (is, "type"); |
|
1489 |
|
1490 if (tag != (char *) NULL && *tag != '\0') |
|
1491 { |
|
1492 char *ptr = strchr (tag, ' '); |
|
1493 if (ptr != (char *) NULL) |
|
1494 { |
|
1495 *ptr = '\0'; |
|
1496 is_global = (strncmp (tag, "global", 6) == 0); |
|
1497 *ptr = ' '; |
|
1498 if (is_global) |
|
1499 ptr++; |
|
1500 else |
|
1501 ptr = tag; |
|
1502 } |
|
1503 else |
|
1504 ptr = tag; |
|
1505 |
|
1506 if (strncmp (ptr, "scalar", 6) == 0) |
|
1507 type_tag = load (is, scalar_constant); |
|
1508 else if (strncmp (ptr, "matrix", 6) == 0) |
|
1509 type_tag = load (is, matrix_constant); |
|
1510 else if (strncmp (ptr, "complex scalar", 14) == 0) |
|
1511 type_tag = load (is, complex_scalar_constant); |
|
1512 else if (strncmp (ptr, "complex matrix", 14) == 0) |
|
1513 type_tag = load (is, complex_matrix_constant); |
|
1514 else if (strncmp (ptr, "string", 6) == 0) |
|
1515 type_tag = load (is, string_constant); |
|
1516 else if (strncmp (ptr, "range", 5) == 0) |
|
1517 type_tag = load (is, range_constant); |
|
1518 else |
|
1519 ::error ("unknown constant type `%s'", tag); |
|
1520 } |
|
1521 else |
|
1522 ::error ("failed to extract keyword specifying value type"); |
|
1523 |
|
1524 delete [] tag; |
|
1525 |
|
1526 return is_global; |
|
1527 } |
|
1528 |
|
1529 tree_constant_rep::constant_type |
|
1530 tree_constant_rep::load (istream& is, tree_constant_rep::constant_type t) |
|
1531 { |
|
1532 tree_constant_rep::constant_type status = unknown_constant; |
|
1533 |
|
1534 switch (t) |
|
1535 { |
|
1536 case scalar_constant: |
|
1537 is >> scalar; |
|
1538 if (is) |
|
1539 status = scalar_constant; |
|
1540 else |
|
1541 ::error ("failed to load scalar constant"); |
|
1542 break; |
|
1543 case matrix_constant: |
|
1544 { |
|
1545 int nr = 0, nc = 0; |
|
1546 |
|
1547 if (extract_keyword (is, "rows", nr) && nr > 0 |
|
1548 && extract_keyword (is, "columns", nc) && nc > 0) |
|
1549 { |
|
1550 matrix = new Matrix (nr, nc); |
|
1551 is >> *matrix; |
|
1552 if (is) |
|
1553 status = matrix_constant; |
|
1554 else |
|
1555 ::error ("failed to load matrix constant"); |
|
1556 } |
|
1557 else |
|
1558 ::error ("failed to extract number of rows and columns"); |
|
1559 } |
|
1560 break; |
|
1561 case complex_scalar_constant: |
|
1562 complex_scalar = new Complex; |
|
1563 is >> *complex_scalar; |
|
1564 if (is) |
|
1565 status = complex_scalar_constant; |
|
1566 else |
|
1567 ::error ("failed to load complex scalar constant"); |
|
1568 break; |
|
1569 case complex_matrix_constant: |
|
1570 { |
|
1571 int nr = 0, nc = 0; |
|
1572 |
|
1573 if (extract_keyword (is, "rows", nr) && nr > 0 |
|
1574 && extract_keyword (is, "columns", nc) && nc > 0) |
|
1575 { |
|
1576 complex_matrix = new ComplexMatrix (nr, nc); |
|
1577 is >> *complex_matrix; |
|
1578 if (is) |
|
1579 status = complex_matrix_constant; |
|
1580 else |
|
1581 ::error ("failed to load complex matrix constant"); |
|
1582 } |
|
1583 else |
|
1584 ::error ("failed to extract number of rows and columns"); |
|
1585 } |
|
1586 break; |
|
1587 case string_constant: |
|
1588 { |
|
1589 int len; |
|
1590 if (extract_keyword (is, "length", len) && len > 0) |
|
1591 { |
|
1592 string = new char [len+1]; |
|
1593 is.get (string, len+1, EOF); |
|
1594 if (is) |
|
1595 status = string_constant; |
|
1596 else |
|
1597 ::error ("failed to load string constant"); |
|
1598 } |
|
1599 else |
|
1600 ::error ("failed to extract string length"); |
|
1601 } |
|
1602 break; |
|
1603 case range_constant: |
|
1604 skip_comments (is); |
|
1605 range = new Range (); |
|
1606 is >> *range; |
|
1607 if (is) |
|
1608 status = range_constant; |
|
1609 else |
|
1610 ::error ("failed to load range constant"); |
|
1611 break; |
|
1612 default: |
|
1613 panic_impossible (); |
|
1614 break; |
|
1615 } |
|
1616 return status; |
|
1617 } |
|
1618 |
|
1619 double |
|
1620 tree_constant_rep::double_value (void) const |
|
1621 { |
|
1622 switch (type_tag) |
|
1623 { |
|
1624 case scalar_constant: |
|
1625 return scalar; |
|
1626 case complex_scalar_constant: |
|
1627 { |
|
1628 int flag = user_pref.ok_to_lose_imaginary_part; |
|
1629 if (flag == -1) |
|
1630 warning ("implicit conversion of complex value to real value"); |
|
1631 |
|
1632 if (flag != 0) |
|
1633 return ::real (*complex_scalar); |
|
1634 |
|
1635 ::error ("implicit conversion of complex value to real value"); |
|
1636 ::error ("not allowed"); |
|
1637 jump_to_top_level (); |
|
1638 } |
|
1639 default: |
|
1640 panic_impossible (); |
|
1641 break; |
|
1642 } |
|
1643 } |
|
1644 |
|
1645 Matrix |
|
1646 tree_constant_rep::matrix_value (void) const |
|
1647 { |
|
1648 switch (type_tag) |
|
1649 { |
|
1650 case scalar_constant: |
|
1651 return Matrix (1, 1, scalar); |
|
1652 case matrix_constant: |
|
1653 return *matrix; |
|
1654 case complex_scalar_constant: |
|
1655 case complex_matrix_constant: |
|
1656 { |
|
1657 int flag = user_pref.ok_to_lose_imaginary_part; |
|
1658 if (flag == -1) |
|
1659 warning ("implicit conversion of complex matrix to real matrix"); |
|
1660 |
|
1661 if (flag != 0) |
|
1662 { |
|
1663 if (type_tag == complex_scalar_constant) |
|
1664 return Matrix (1, 1, ::real (*complex_scalar)); |
|
1665 else if (type_tag == complex_matrix_constant) |
|
1666 return ::real (*complex_matrix); |
|
1667 else |
|
1668 panic_impossible (); |
|
1669 } |
|
1670 else |
|
1671 { |
|
1672 ::error ("implicit conversion of complex matrix to real matrix"); |
|
1673 ::error ("not allowed"); |
|
1674 } |
|
1675 jump_to_top_level (); |
|
1676 } |
|
1677 default: |
|
1678 panic_impossible (); |
|
1679 break; |
|
1680 } |
|
1681 } |
|
1682 |
|
1683 Complex |
|
1684 tree_constant_rep::complex_value (void) const |
|
1685 { |
|
1686 switch (type_tag) |
|
1687 { |
|
1688 case complex_scalar_constant: |
|
1689 return *complex_scalar; |
|
1690 case scalar_constant: |
|
1691 return Complex (scalar); |
|
1692 default: |
|
1693 panic_impossible (); |
|
1694 break; |
|
1695 } |
|
1696 } |
|
1697 |
|
1698 ComplexMatrix |
|
1699 tree_constant_rep::complex_matrix_value (void) const |
|
1700 { |
|
1701 switch (type_tag) |
|
1702 { |
|
1703 case scalar_constant: |
|
1704 { |
|
1705 return ComplexMatrix (1, 1, Complex (scalar)); |
|
1706 } |
|
1707 case complex_scalar_constant: |
|
1708 { |
|
1709 return ComplexMatrix (1, 1, *complex_scalar); |
|
1710 } |
|
1711 case matrix_constant: |
|
1712 { |
|
1713 return ComplexMatrix (*matrix); |
|
1714 } |
|
1715 case complex_matrix_constant: |
|
1716 return *complex_matrix; |
|
1717 break; |
|
1718 default: |
|
1719 panic_impossible (); |
|
1720 break; |
|
1721 } |
|
1722 } |
|
1723 |
|
1724 char * |
|
1725 tree_constant_rep::string_value (void) const |
|
1726 { |
|
1727 assert (type_tag == string_constant); |
|
1728 return string; |
|
1729 } |
|
1730 |
|
1731 Range |
|
1732 tree_constant_rep::range_value (void) const |
|
1733 { |
|
1734 assert (type_tag == range_constant); |
|
1735 return *range; |
|
1736 } |
|
1737 |
|
1738 int |
|
1739 tree_constant_rep::rows (void) const |
|
1740 { |
|
1741 int retval = -1; |
|
1742 switch (type_tag) |
|
1743 { |
|
1744 case scalar_constant: |
|
1745 case complex_scalar_constant: |
|
1746 retval = 1; |
|
1747 break; |
|
1748 case string_constant: |
|
1749 case range_constant: |
|
1750 retval = (columns () > 0); |
|
1751 break; |
|
1752 case matrix_constant: |
|
1753 retval = matrix->rows (); |
|
1754 break; |
|
1755 case complex_matrix_constant: |
|
1756 retval = complex_matrix->rows (); |
|
1757 break; |
|
1758 case magic_colon: |
|
1759 ::error ("invalid use of colon operator"); |
|
1760 break; |
|
1761 case unknown_constant: |
|
1762 retval = 0; |
|
1763 break; |
|
1764 default: |
|
1765 panic_impossible (); |
|
1766 break; |
|
1767 } |
|
1768 return retval; |
|
1769 } |
|
1770 |
|
1771 int |
|
1772 tree_constant_rep::columns (void) const |
|
1773 { |
|
1774 int retval = -1; |
|
1775 switch (type_tag) |
|
1776 { |
|
1777 case scalar_constant: |
|
1778 case complex_scalar_constant: |
|
1779 retval = 1; |
|
1780 break; |
|
1781 case matrix_constant: |
|
1782 retval = matrix->columns (); |
|
1783 break; |
|
1784 case complex_matrix_constant: |
|
1785 retval = complex_matrix->columns (); |
|
1786 break; |
|
1787 case string_constant: |
|
1788 retval = strlen (string); |
|
1789 break; |
|
1790 case range_constant: |
|
1791 retval = range->nelem (); |
|
1792 break; |
|
1793 case magic_colon: |
|
1794 ::error ("invalid use of colon operator"); |
|
1795 break; |
|
1796 case unknown_constant: |
|
1797 retval = 0; |
|
1798 break; |
|
1799 default: |
|
1800 panic_impossible (); |
|
1801 break; |
|
1802 } |
|
1803 return retval; |
|
1804 } |
|
1805 |
|
1806 tree_constant |
|
1807 tree_constant_rep::all (void) const |
|
1808 { |
|
1809 if (type_tag == string_constant || type_tag == range_constant) |
|
1810 { |
|
1811 tree_constant tmp = make_numeric (); |
|
1812 return tmp.all (); |
|
1813 } |
|
1814 |
|
1815 tree_constant retval; |
|
1816 switch (type_tag) |
|
1817 { |
|
1818 case scalar_constant: |
|
1819 { |
|
1820 double status = (scalar != 0.0); |
|
1821 retval = tree_constant (status); |
|
1822 } |
|
1823 break; |
|
1824 case matrix_constant: |
|
1825 { |
|
1826 Matrix m = matrix->all (); |
|
1827 retval = tree_constant (m); |
|
1828 } |
|
1829 break; |
|
1830 case complex_scalar_constant: |
|
1831 { |
|
1832 double status = (*complex_scalar != 0.0); |
|
1833 retval = tree_constant (status); |
|
1834 } |
|
1835 break; |
|
1836 case complex_matrix_constant: |
|
1837 { |
|
1838 Matrix m = complex_matrix->all (); |
|
1839 retval = tree_constant (m); |
|
1840 } |
|
1841 break; |
|
1842 case string_constant: |
|
1843 case range_constant: |
|
1844 case magic_colon: |
|
1845 default: |
|
1846 panic_impossible (); |
|
1847 break; |
|
1848 } |
|
1849 return retval; |
|
1850 } |
|
1851 |
|
1852 tree_constant |
|
1853 tree_constant_rep::any (void) const |
|
1854 { |
|
1855 if (type_tag == string_constant || type_tag == range_constant) |
|
1856 { |
|
1857 tree_constant tmp = make_numeric (); |
|
1858 return tmp.any (); |
|
1859 } |
|
1860 |
|
1861 tree_constant retval; |
|
1862 switch (type_tag) |
|
1863 { |
|
1864 case scalar_constant: |
|
1865 { |
|
1866 double status = (scalar != 0.0); |
|
1867 retval = tree_constant (status); |
|
1868 } |
|
1869 break; |
|
1870 case matrix_constant: |
|
1871 { |
|
1872 Matrix m = matrix->any (); |
|
1873 retval = tree_constant (m); |
|
1874 } |
|
1875 break; |
|
1876 case complex_scalar_constant: |
|
1877 { |
|
1878 double status = (*complex_scalar != 0.0); |
|
1879 retval = tree_constant (status); |
|
1880 } |
|
1881 break; |
|
1882 case complex_matrix_constant: |
|
1883 { |
|
1884 Matrix m = complex_matrix->any (); |
|
1885 retval = tree_constant (m); |
|
1886 } |
|
1887 break; |
|
1888 case string_constant: |
|
1889 case range_constant: |
|
1890 case magic_colon: |
|
1891 default: |
|
1892 panic_impossible (); |
|
1893 break; |
|
1894 } |
|
1895 return retval; |
|
1896 } |
|
1897 |
|
1898 tree_constant |
|
1899 tree_constant_rep::isstr (void) const |
|
1900 { |
|
1901 double status = 0.0; |
|
1902 if (type_tag == string_constant) |
|
1903 status = 1.0; |
|
1904 tree_constant retval (status); |
|
1905 return retval; |
|
1906 } |
|
1907 |
|
1908 tree_constant |
|
1909 tree_constant_rep::convert_to_str (void) |
|
1910 { |
|
1911 tree_constant retval; |
|
1912 |
|
1913 switch (type_tag) |
|
1914 { |
|
1915 case complex_scalar_constant: |
|
1916 case scalar_constant: |
|
1917 { |
|
1918 double d = double_value (); |
|
1919 int i = NINT (d); |
|
1920 // Warn about out of range conversions? |
|
1921 char s[2]; |
|
1922 s[0] = (char) i; |
|
1923 s[1] = '\0'; |
|
1924 retval = tree_constant (s); |
|
1925 } |
|
1926 break; |
|
1927 case complex_matrix_constant: |
|
1928 case matrix_constant: |
|
1929 { |
|
1930 ColumnVector v = to_vector (); |
|
1931 int len = v.length (); |
|
1932 if (len == 0) |
|
1933 ::error ("can only convert vectors and scalars to strings"); |
|
1934 else |
|
1935 { |
|
1936 char *s = new char [len+1]; |
|
1937 s[len] = '\0'; |
|
1938 for (int i = 0; i < len; i++) |
|
1939 { |
|
1940 double d = v.elem (i); |
|
1941 int ival = NINT (d); |
|
1942 // Warn about out of range conversions? |
|
1943 s[i] = (char) ival; |
|
1944 } |
|
1945 retval = tree_constant (s); |
|
1946 delete [] s; |
|
1947 } |
|
1948 } |
|
1949 break; |
|
1950 case range_constant: |
|
1951 { |
|
1952 Range r = range_value (); |
|
1953 double b = r.base (); |
|
1954 double incr = r.inc (); |
|
1955 int nel = r.nelem (); |
|
1956 char *s = new char [nel+1]; |
|
1957 s[nel] = '\0'; |
|
1958 for (int i = 0; i < nel; i++) |
|
1959 { |
|
1960 double d = b + i * incr; |
|
1961 int ival = NINT (d); |
|
1962 // Warn about out of range conversions? |
|
1963 s[i] = (char) ival; |
|
1964 } |
|
1965 retval = tree_constant (s); |
|
1966 delete [] s; |
|
1967 } |
|
1968 break; |
|
1969 case string_constant: |
|
1970 retval = tree_constant (*this); |
|
1971 break; |
|
1972 case magic_colon: |
|
1973 default: |
|
1974 panic_impossible (); |
|
1975 break; |
|
1976 } |
|
1977 return retval; |
|
1978 } |
|
1979 |
|
1980 void |
|
1981 tree_constant_rep::convert_to_row_or_column_vector (void) |
|
1982 { |
|
1983 assert (type_tag == matrix_constant || type_tag == complex_matrix_constant); |
|
1984 |
|
1985 int nr = rows (); |
|
1986 int nc = columns (); |
|
1987 |
|
1988 int len = nr * nc; |
|
1989 |
|
1990 assert (len > 0); |
|
1991 |
|
1992 int new_nr = 1; |
|
1993 int new_nc = 1; |
|
1994 |
|
1995 if (user_pref.prefer_column_vectors) |
|
1996 new_nr = len; |
|
1997 else |
|
1998 new_nc = len; |
|
1999 |
|
2000 if (type_tag == matrix_constant) |
|
2001 { |
|
2002 Matrix *m = new Matrix (new_nr, new_nc); |
|
2003 |
|
2004 double *cop_out = matrix->fortran_vec (); |
|
2005 |
|
2006 for (int i = 0; i < len; i++) |
|
2007 { |
|
2008 if (new_nr == 1) |
|
2009 m->elem (0, i) = *cop_out++; |
|
2010 else |
|
2011 m->elem (i, 0) = *cop_out++; |
|
2012 } |
|
2013 |
|
2014 delete matrix; |
|
2015 matrix = m; |
|
2016 } |
|
2017 else |
|
2018 { |
|
2019 ComplexMatrix *cm = new ComplexMatrix (new_nr, new_nc); |
|
2020 |
|
2021 Complex *cop_out = complex_matrix->fortran_vec (); |
|
2022 |
|
2023 for (int i = 0; i < len; i++) |
|
2024 { |
|
2025 if (new_nr == 1) |
|
2026 cm->elem (0, i) = *cop_out++; |
|
2027 else |
|
2028 cm->elem (i, 0) = *cop_out++; |
|
2029 } |
|
2030 |
|
2031 delete complex_matrix; |
|
2032 complex_matrix = cm; |
|
2033 } |
|
2034 } |
|
2035 |
|
2036 int |
|
2037 tree_constant_rep::is_true (void) const |
|
2038 { |
|
2039 if (type_tag == string_constant || type_tag == range_constant) |
|
2040 { |
|
2041 tree_constant tmp = make_numeric (); |
|
2042 return tmp.is_true (); |
|
2043 } |
|
2044 |
|
2045 int retval; |
|
2046 switch (type_tag) |
|
2047 { |
|
2048 case scalar_constant: |
|
2049 retval = (scalar != 0.0); |
|
2050 break; |
|
2051 case matrix_constant: |
|
2052 { |
|
2053 Matrix m = (matrix->all ()) . all (); |
|
2054 retval = (m.rows () == 1 |
|
2055 && m.columns () == 1 |
|
2056 && m.elem (0, 0) != 0.0); |
|
2057 } |
|
2058 break; |
|
2059 case complex_scalar_constant: |
|
2060 retval = (*complex_scalar != 0.0); |
|
2061 break; |
|
2062 case complex_matrix_constant: |
|
2063 { |
|
2064 Matrix m = (complex_matrix->all ()) . all (); |
|
2065 retval = (m.rows () == 1 |
|
2066 && m.columns () == 1 |
|
2067 && m.elem (0, 0) != 0.0); |
|
2068 } |
|
2069 break; |
|
2070 case string_constant: |
|
2071 case range_constant: |
|
2072 case magic_colon: |
|
2073 default: |
|
2074 panic_impossible (); |
|
2075 break; |
|
2076 } |
|
2077 return retval; |
|
2078 } |
|
2079 |
|
2080 tree_constant |
|
2081 tree_constant_rep::cumprod (void) const |
|
2082 { |
|
2083 if (type_tag == string_constant || type_tag == range_constant) |
|
2084 { |
|
2085 tree_constant tmp = make_numeric (); |
|
2086 return tmp.cumprod (); |
|
2087 } |
|
2088 |
|
2089 tree_constant retval; |
|
2090 switch (type_tag) |
|
2091 { |
|
2092 case scalar_constant: |
|
2093 retval = tree_constant (scalar); |
|
2094 break; |
|
2095 case matrix_constant: |
|
2096 { |
|
2097 Matrix m = matrix->cumprod (); |
|
2098 retval = tree_constant (m); |
|
2099 } |
|
2100 break; |
|
2101 case complex_scalar_constant: |
|
2102 retval = tree_constant (*complex_scalar); |
|
2103 break; |
|
2104 case complex_matrix_constant: |
|
2105 { |
|
2106 ComplexMatrix m = complex_matrix->cumprod (); |
|
2107 retval = tree_constant (m); |
|
2108 } |
|
2109 break; |
|
2110 case string_constant: |
|
2111 case range_constant: |
|
2112 case magic_colon: |
|
2113 default: |
|
2114 panic_impossible (); |
|
2115 break; |
|
2116 } |
|
2117 return retval; |
|
2118 } |
|
2119 |
|
2120 tree_constant |
|
2121 tree_constant_rep::cumsum (void) const |
|
2122 { |
|
2123 if (type_tag == string_constant || type_tag == range_constant) |
|
2124 { |
|
2125 tree_constant tmp = make_numeric (); |
|
2126 return tmp.cumsum (); |
|
2127 } |
|
2128 |
|
2129 tree_constant retval; |
|
2130 switch (type_tag) |
|
2131 { |
|
2132 case scalar_constant: |
|
2133 retval = tree_constant (scalar); |
|
2134 break; |
|
2135 case matrix_constant: |
|
2136 { |
|
2137 Matrix m = matrix->cumsum (); |
|
2138 retval = tree_constant (m); |
|
2139 } |
|
2140 break; |
|
2141 case complex_scalar_constant: |
|
2142 retval = tree_constant (*complex_scalar); |
|
2143 break; |
|
2144 case complex_matrix_constant: |
|
2145 { |
|
2146 ComplexMatrix m = complex_matrix->cumsum (); |
|
2147 retval = tree_constant (m); |
|
2148 } |
|
2149 break; |
|
2150 case string_constant: |
|
2151 case range_constant: |
|
2152 case magic_colon: |
|
2153 default: |
|
2154 panic_impossible (); |
|
2155 break; |
|
2156 } |
|
2157 return retval; |
|
2158 } |
|
2159 |
|
2160 tree_constant |
|
2161 tree_constant_rep::prod (void) const |
|
2162 { |
|
2163 if (type_tag == string_constant || type_tag == range_constant) |
|
2164 { |
|
2165 tree_constant tmp = make_numeric (); |
|
2166 return tmp.prod (); |
|
2167 } |
|
2168 |
|
2169 tree_constant retval; |
|
2170 switch (type_tag) |
|
2171 { |
|
2172 case scalar_constant: |
|
2173 retval = tree_constant (scalar); |
|
2174 break; |
|
2175 case matrix_constant: |
|
2176 { |
|
2177 Matrix m = matrix->prod (); |
|
2178 retval = tree_constant (m); |
|
2179 } |
|
2180 break; |
|
2181 case complex_scalar_constant: |
|
2182 retval = tree_constant (*complex_scalar); |
|
2183 break; |
|
2184 case complex_matrix_constant: |
|
2185 { |
|
2186 ComplexMatrix m = complex_matrix->prod (); |
|
2187 retval = tree_constant (m); |
|
2188 } |
|
2189 break; |
|
2190 case string_constant: |
|
2191 case range_constant: |
|
2192 case magic_colon: |
|
2193 default: |
|
2194 panic_impossible (); |
|
2195 break; |
|
2196 } |
|
2197 return retval; |
|
2198 } |
|
2199 |
|
2200 tree_constant |
|
2201 tree_constant_rep::sum (void) const |
|
2202 { |
|
2203 if (type_tag == string_constant || type_tag == range_constant) |
|
2204 { |
|
2205 tree_constant tmp = make_numeric (); |
|
2206 return tmp.sum (); |
|
2207 } |
|
2208 |
|
2209 tree_constant retval; |
|
2210 switch (type_tag) |
|
2211 { |
|
2212 case scalar_constant: |
|
2213 retval = tree_constant (scalar); |
|
2214 break; |
|
2215 case matrix_constant: |
|
2216 { |
|
2217 Matrix m = matrix->sum (); |
|
2218 retval = tree_constant (m); |
|
2219 } |
|
2220 break; |
|
2221 case complex_scalar_constant: |
|
2222 retval = tree_constant (*complex_scalar); |
|
2223 break; |
|
2224 case complex_matrix_constant: |
|
2225 { |
|
2226 ComplexMatrix m = complex_matrix->sum (); |
|
2227 retval = tree_constant (m); |
|
2228 } |
|
2229 break; |
|
2230 case string_constant: |
|
2231 case range_constant: |
|
2232 case magic_colon: |
|
2233 default: |
|
2234 panic_impossible (); |
|
2235 break; |
|
2236 } |
|
2237 return retval; |
|
2238 } |
|
2239 |
|
2240 tree_constant |
|
2241 tree_constant_rep::sumsq (void) const |
|
2242 { |
|
2243 if (type_tag == string_constant || type_tag == range_constant) |
|
2244 { |
|
2245 tree_constant tmp = make_numeric (); |
|
2246 return tmp.sumsq (); |
|
2247 } |
|
2248 |
|
2249 tree_constant retval; |
|
2250 switch (type_tag) |
|
2251 { |
|
2252 case scalar_constant: |
|
2253 retval = tree_constant (scalar * scalar); |
|
2254 break; |
|
2255 case matrix_constant: |
|
2256 { |
|
2257 Matrix m = matrix->sumsq (); |
|
2258 retval = tree_constant (m); |
|
2259 } |
|
2260 break; |
|
2261 case complex_scalar_constant: |
|
2262 { |
|
2263 Complex c (*complex_scalar); |
|
2264 retval = tree_constant (c * c); |
|
2265 } |
|
2266 break; |
|
2267 case complex_matrix_constant: |
|
2268 { |
|
2269 ComplexMatrix m = complex_matrix->sumsq (); |
|
2270 retval = tree_constant (m); |
|
2271 } |
|
2272 break; |
|
2273 case string_constant: |
|
2274 case range_constant: |
|
2275 case magic_colon: |
|
2276 default: |
|
2277 panic_impossible (); |
|
2278 break; |
|
2279 } |
|
2280 return retval; |
|
2281 } |
|
2282 |
|
2283 static tree_constant |
|
2284 make_diag (const Matrix& v, int k) |
|
2285 { |
|
2286 int nr = v.rows (); |
|
2287 int nc = v.columns (); |
|
2288 assert (nc == 1 || nr == 1); |
|
2289 |
|
2290 tree_constant retval; |
|
2291 |
|
2292 int roff = 0; |
|
2293 int coff = 0; |
|
2294 if (k > 0) |
|
2295 { |
|
2296 roff = 0; |
|
2297 coff = k; |
|
2298 } |
|
2299 else if (k < 0) |
|
2300 { |
|
2301 roff = -k; |
|
2302 coff = 0; |
|
2303 } |
|
2304 |
|
2305 if (nr == 1) |
|
2306 { |
|
2307 int n = nc + ABS (k); |
|
2308 Matrix m (n, n, 0.0); |
|
2309 for (int i = 0; i < nc; i++) |
|
2310 m.elem (i+roff, i+coff) = v.elem (0, i); |
|
2311 retval = tree_constant (m); |
|
2312 } |
|
2313 else |
|
2314 { |
|
2315 int n = nr + ABS (k); |
|
2316 Matrix m (n, n, 0.0); |
|
2317 for (int i = 0; i < nr; i++) |
|
2318 m.elem (i+roff, i+coff) = v.elem (i, 0); |
|
2319 retval = tree_constant (m); |
|
2320 } |
|
2321 |
|
2322 return retval; |
|
2323 } |
|
2324 |
|
2325 static tree_constant |
|
2326 make_diag (const ComplexMatrix& v, int k) |
|
2327 { |
|
2328 int nr = v.rows (); |
|
2329 int nc = v.columns (); |
|
2330 assert (nc == 1 || nr == 1); |
|
2331 |
|
2332 tree_constant retval; |
|
2333 |
|
2334 int roff = 0; |
|
2335 int coff = 0; |
|
2336 if (k > 0) |
|
2337 { |
|
2338 roff = 0; |
|
2339 coff = k; |
|
2340 } |
|
2341 else if (k < 0) |
|
2342 { |
|
2343 roff = -k; |
|
2344 coff = 0; |
|
2345 } |
|
2346 |
|
2347 if (nr == 1) |
|
2348 { |
|
2349 int n = nc + ABS (k); |
|
2350 ComplexMatrix m (n, n, 0.0); |
|
2351 for (int i = 0; i < nc; i++) |
|
2352 m.elem (i+roff, i+coff) = v.elem (0, i); |
|
2353 retval = tree_constant (m); |
|
2354 } |
|
2355 else |
|
2356 { |
|
2357 int n = nr + ABS (k); |
|
2358 ComplexMatrix m (n, n, 0.0); |
|
2359 for (int i = 0; i < nr; i++) |
|
2360 m.elem (i+roff, i+coff) = v.elem (i, 0); |
|
2361 retval = tree_constant (m); |
|
2362 } |
|
2363 |
|
2364 return retval; |
|
2365 } |
|
2366 |
|
2367 tree_constant |
|
2368 tree_constant_rep::diag (void) const |
|
2369 { |
|
2370 if (type_tag == string_constant || type_tag == range_constant) |
|
2371 { |
|
2372 tree_constant tmp = make_numeric (); |
|
2373 return tmp.diag (); |
|
2374 } |
|
2375 |
|
2376 tree_constant retval; |
|
2377 switch (type_tag) |
|
2378 { |
|
2379 case scalar_constant: |
|
2380 retval = tree_constant (scalar); |
|
2381 break; |
|
2382 case matrix_constant: |
|
2383 { |
|
2384 int nr = rows (); |
|
2385 int nc = columns (); |
|
2386 if (nr == 0 || nc == 0) |
|
2387 { |
|
2388 Matrix mtmp; |
|
2389 retval = tree_constant (mtmp); |
|
2390 } |
|
2391 else if (nr == 1 || nc == 1) |
|
2392 retval = make_diag (matrix_value (), 0); |
|
2393 else |
|
2394 { |
|
2395 ColumnVector v = matrix->diag (); |
|
2396 if (v.capacity () > 0) |
|
2397 retval = tree_constant (v); |
|
2398 } |
|
2399 } |
|
2400 break; |
|
2401 case complex_scalar_constant: |
|
2402 retval = tree_constant (*complex_scalar); |
|
2403 break; |
|
2404 case complex_matrix_constant: |
|
2405 { |
|
2406 int nr = rows (); |
|
2407 int nc = columns (); |
|
2408 if (nr == 0 || nc == 0) |
|
2409 { |
|
2410 Matrix mtmp; |
|
2411 retval = tree_constant (mtmp); |
|
2412 } |
|
2413 else if (nr == 1 || nc == 1) |
|
2414 retval = make_diag (complex_matrix_value (), 0); |
|
2415 else |
|
2416 { |
|
2417 ComplexColumnVector v = complex_matrix->diag (); |
|
2418 if (v.capacity () > 0) |
|
2419 retval = tree_constant (v); |
|
2420 } |
|
2421 } |
|
2422 break; |
|
2423 case string_constant: |
|
2424 case range_constant: |
|
2425 case magic_colon: |
|
2426 default: |
|
2427 panic_impossible (); |
|
2428 break; |
|
2429 } |
|
2430 return retval; |
|
2431 } |
|
2432 |
|
2433 tree_constant |
|
2434 tree_constant_rep::diag (const tree_constant& a) const |
|
2435 { |
|
2436 if (type_tag == string_constant || type_tag == range_constant) |
|
2437 { |
|
2438 tree_constant tmp = make_numeric (); |
|
2439 return tmp.diag (a); |
|
2440 } |
|
2441 |
|
2442 tree_constant tmp_a = a.make_numeric (); |
|
2443 |
|
2444 tree_constant_rep::constant_type a_type = tmp_a.const_type (); |
|
2445 |
|
2446 tree_constant retval; |
|
2447 |
|
2448 switch (type_tag) |
|
2449 { |
|
2450 case scalar_constant: |
|
2451 if (a_type == scalar_constant) |
|
2452 { |
|
2453 int k = NINT (tmp_a.double_value ()); |
|
2454 int n = ABS (k) + 1; |
|
2455 if (k == 0) |
|
2456 retval = tree_constant (scalar); |
|
2457 else if (k > 0) |
|
2458 { |
|
2459 Matrix m (n, n, 0.0); |
|
2460 m.elem (0, k) = scalar; |
|
2461 retval = tree_constant (m); |
|
2462 } |
|
2463 else if (k < 0) |
|
2464 { |
|
2465 Matrix m (n, n, 0.0); |
|
2466 m.elem (-k, 0) = scalar; |
|
2467 retval = tree_constant (m); |
|
2468 } |
|
2469 } |
|
2470 break; |
|
2471 case matrix_constant: |
|
2472 if (a_type == scalar_constant) |
|
2473 { |
|
2474 int k = NINT (tmp_a.double_value ()); |
|
2475 int nr = rows (); |
|
2476 int nc = columns (); |
|
2477 if (nr == 0 || nc == 0) |
|
2478 { |
|
2479 Matrix mtmp; |
|
2480 retval = tree_constant (mtmp); |
|
2481 } |
|
2482 else if (nr == 1 || nc == 1) |
|
2483 retval = make_diag (matrix_value (), k); |
|
2484 else |
|
2485 { |
|
2486 ColumnVector d = matrix->diag (k); |
|
2487 retval = tree_constant (d); |
|
2488 } |
|
2489 } |
|
2490 else |
|
2491 ::error ("diag: invalid second argument"); |
|
2492 |
|
2493 break; |
|
2494 case complex_scalar_constant: |
|
2495 if (a_type == scalar_constant) |
|
2496 { |
|
2497 int k = NINT (tmp_a.double_value ()); |
|
2498 int n = ABS (k) + 1; |
|
2499 if (k == 0) |
|
2500 retval = tree_constant (*complex_scalar); |
|
2501 else if (k > 0) |
|
2502 { |
|
2503 ComplexMatrix m (n, n, 0.0); |
|
2504 m.elem (0, k) = *complex_scalar; |
|
2505 retval = tree_constant (m); |
|
2506 } |
|
2507 else if (k < 0) |
|
2508 { |
|
2509 ComplexMatrix m (n, n, 0.0); |
|
2510 m.elem (-k, 0) = *complex_scalar; |
|
2511 retval = tree_constant (m); |
|
2512 } |
|
2513 } |
|
2514 break; |
|
2515 case complex_matrix_constant: |
|
2516 if (a_type == scalar_constant) |
|
2517 { |
|
2518 int k = NINT (tmp_a.double_value ()); |
|
2519 int nr = rows (); |
|
2520 int nc = columns (); |
|
2521 if (nr == 0 || nc == 0) |
|
2522 { |
|
2523 Matrix mtmp; |
|
2524 retval = tree_constant (mtmp); |
|
2525 } |
|
2526 else if (nr == 1 || nc == 1) |
|
2527 retval = make_diag (complex_matrix_value (), k); |
|
2528 else |
|
2529 { |
|
2530 ComplexColumnVector d = complex_matrix->diag (k); |
|
2531 retval = tree_constant (d); |
|
2532 } |
|
2533 } |
|
2534 else |
|
2535 ::error ("diag: invalid second argument"); |
|
2536 |
|
2537 break; |
|
2538 case string_constant: |
|
2539 case range_constant: |
|
2540 case magic_colon: |
|
2541 default: |
|
2542 panic_impossible (); |
|
2543 break; |
|
2544 } |
|
2545 return retval; |
|
2546 } |
|
2547 |
|
2548 tree_constant |
|
2549 tree_constant_rep::mapper (Mapper_fcn& m_fcn, int print) const |
|
2550 { |
|
2551 tree_constant retval; |
|
2552 |
|
2553 if (type_tag == string_constant || type_tag == range_constant) |
|
2554 { |
|
2555 tree_constant tmp = make_numeric (); |
|
2556 return tmp.mapper (m_fcn, print); |
|
2557 } |
|
2558 |
|
2559 switch (type_tag) |
|
2560 { |
|
2561 case scalar_constant: |
|
2562 if (m_fcn.can_return_complex_for_real_arg |
|
2563 && (scalar < m_fcn.lower_limit |
|
2564 || scalar > m_fcn.upper_limit)) |
|
2565 { |
|
2566 if (m_fcn.c_c_mapper != NULL) |
|
2567 { |
|
2568 Complex c = m_fcn.c_c_mapper (Complex (scalar)); |
|
2569 retval = tree_constant (c); |
|
2570 } |
|
2571 else |
|
2572 panic_impossible (); |
|
2573 } |
|
2574 else |
|
2575 { |
|
2576 if (m_fcn.d_d_mapper != NULL) |
|
2577 { |
|
2578 double d = m_fcn.d_d_mapper (scalar); |
|
2579 retval = tree_constant (d); |
|
2580 } |
|
2581 else |
|
2582 panic_impossible (); |
|
2583 } |
|
2584 break; |
|
2585 case matrix_constant: |
|
2586 if (m_fcn.can_return_complex_for_real_arg |
|
2587 && (any_element_less_than (*matrix, m_fcn.lower_limit) |
|
2588 || any_element_greater_than (*matrix, m_fcn.upper_limit))) |
|
2589 { |
|
2590 if (m_fcn.c_c_mapper != NULL) |
|
2591 { |
|
2592 ComplexMatrix cm = map (m_fcn.c_c_mapper, |
|
2593 ComplexMatrix (*matrix)); |
|
2594 retval = tree_constant (cm); |
|
2595 } |
|
2596 else |
|
2597 panic_impossible (); |
|
2598 } |
|
2599 else |
|
2600 { |
|
2601 if (m_fcn.d_d_mapper != NULL) |
|
2602 { |
|
2603 Matrix m = map (m_fcn.d_d_mapper, *matrix); |
|
2604 retval = tree_constant (m); |
|
2605 } |
|
2606 else |
|
2607 panic_impossible (); |
|
2608 } |
|
2609 break; |
|
2610 case complex_scalar_constant: |
|
2611 if (m_fcn.d_c_mapper != NULL) |
|
2612 { |
|
2613 double d; |
|
2614 d = m_fcn.d_c_mapper (*complex_scalar); |
|
2615 retval = tree_constant (d); |
|
2616 } |
|
2617 else if (m_fcn.c_c_mapper != NULL) |
|
2618 { |
|
2619 Complex c; |
|
2620 c = m_fcn.c_c_mapper (*complex_scalar); |
|
2621 retval = tree_constant (c); |
|
2622 } |
|
2623 else |
|
2624 panic_impossible (); |
|
2625 break; |
|
2626 case complex_matrix_constant: |
|
2627 if (m_fcn.d_c_mapper != NULL) |
|
2628 { |
|
2629 Matrix m; |
|
2630 m = map (m_fcn.d_c_mapper, *complex_matrix); |
|
2631 retval = tree_constant (m); |
|
2632 } |
|
2633 else if (m_fcn.c_c_mapper != NULL) |
|
2634 { |
|
2635 ComplexMatrix cm; |
|
2636 cm = map (m_fcn.c_c_mapper, *complex_matrix); |
|
2637 retval = tree_constant (cm); |
|
2638 } |
|
2639 else |
|
2640 panic_impossible (); |
|
2641 break; |
|
2642 case string_constant: |
|
2643 case range_constant: |
|
2644 case magic_colon: |
|
2645 default: |
|
2646 panic_impossible (); |
|
2647 break; |
|
2648 } |
|
2649 return retval; |
|
2650 } |
|
2651 |
|
2652 /* |
|
2653 * Top-level tree-constant function that handles assignments. Only |
|
2654 * decide if the left-hand side is currently a scalar or a matrix and |
|
2655 * hand off to other functions to do the real work. |
|
2656 */ |
|
2657 void |
500
|
2658 tree_constant_rep::assign (tree_constant& rhs, |
|
2659 const Octave_object& args, int nargs) |
492
|
2660 { |
|
2661 tree_constant rhs_tmp = rhs.make_numeric (); |
|
2662 |
|
2663 // This is easier than actually handling assignments to strings. |
|
2664 // An assignment to a range will normally require a conversion to a |
|
2665 // vector since it will normally destroy the equally-spaced property |
|
2666 // of the range elements. |
|
2667 |
|
2668 if (type_tag == string_constant || type_tag == range_constant) |
|
2669 force_numeric (); |
|
2670 |
|
2671 switch (type_tag) |
|
2672 { |
|
2673 case complex_scalar_constant: |
|
2674 case scalar_constant: |
|
2675 case unknown_constant: |
|
2676 do_scalar_assignment (rhs_tmp, args, nargs); |
|
2677 break; |
|
2678 case complex_matrix_constant: |
|
2679 case matrix_constant: |
|
2680 do_matrix_assignment (rhs_tmp, args, nargs); |
|
2681 break; |
|
2682 case string_constant: |
|
2683 ::error ("invalid assignment to string type"); |
|
2684 break; |
|
2685 case range_constant: |
|
2686 case magic_colon: |
|
2687 default: |
|
2688 panic_impossible (); |
|
2689 break; |
|
2690 } |
|
2691 } |
|
2692 |
|
2693 /* |
|
2694 * Assignments to scalars. If resize_on_range_error is true, |
|
2695 * this can convert the left-hand side to a matrix. |
|
2696 */ |
|
2697 void |
|
2698 tree_constant_rep::do_scalar_assignment (tree_constant& rhs, |
500
|
2699 const Octave_object& args, int nargs) |
492
|
2700 { |
|
2701 assert (type_tag == unknown_constant |
|
2702 || type_tag == scalar_constant |
|
2703 || type_tag == complex_scalar_constant); |
|
2704 |
|
2705 if ((rhs.is_scalar_type () || rhs.is_zero_by_zero ()) |
|
2706 && valid_scalar_indices (args, nargs)) |
|
2707 { |
|
2708 if (rhs.is_zero_by_zero ()) |
|
2709 { |
|
2710 if (type_tag == complex_scalar_constant) |
|
2711 delete complex_scalar; |
|
2712 |
|
2713 matrix = new Matrix (0, 0); |
|
2714 type_tag = matrix_constant; |
|
2715 } |
|
2716 else if (type_tag == unknown_constant || type_tag == scalar_constant) |
|
2717 { |
|
2718 if (rhs.const_type () == scalar_constant) |
|
2719 { |
|
2720 scalar = rhs.double_value (); |
|
2721 type_tag = scalar_constant; |
|
2722 } |
|
2723 else if (rhs.const_type () == complex_scalar_constant) |
|
2724 { |
|
2725 complex_scalar = new Complex (rhs.complex_value ()); |
|
2726 type_tag = complex_scalar_constant; |
|
2727 } |
|
2728 else |
|
2729 { |
|
2730 ::error ("invalid assignment to scalar"); |
|
2731 return; |
|
2732 } |
|
2733 } |
|
2734 else |
|
2735 { |
|
2736 if (rhs.const_type () == scalar_constant) |
|
2737 { |
|
2738 delete complex_scalar; |
|
2739 scalar = rhs.double_value (); |
|
2740 type_tag = scalar_constant; |
|
2741 } |
|
2742 else if (rhs.const_type () == complex_scalar_constant) |
|
2743 { |
|
2744 *complex_scalar = rhs.complex_value (); |
|
2745 type_tag = complex_scalar_constant; |
|
2746 } |
|
2747 else |
|
2748 { |
|
2749 ::error ("invalid assignment to scalar"); |
|
2750 return; |
|
2751 } |
|
2752 } |
|
2753 } |
|
2754 else if (user_pref.resize_on_range_error) |
|
2755 { |
|
2756 tree_constant_rep::constant_type old_type_tag = type_tag; |
|
2757 |
|
2758 if (type_tag == complex_scalar_constant) |
|
2759 { |
|
2760 Complex *old_complex = complex_scalar; |
|
2761 complex_matrix = new ComplexMatrix (1, 1, *complex_scalar); |
|
2762 type_tag = complex_matrix_constant; |
|
2763 delete old_complex; |
|
2764 } |
|
2765 else if (type_tag == scalar_constant) |
|
2766 { |
|
2767 matrix = new Matrix (1, 1, scalar); |
|
2768 type_tag = matrix_constant; |
|
2769 } |
|
2770 |
|
2771 // If there is an error, the call to do_matrix_assignment should not |
|
2772 // destroy the current value. tree_constant_rep::eval(int) will take |
|
2773 // care of converting single element matrices back to scalars. |
|
2774 |
|
2775 do_matrix_assignment (rhs, args, nargs); |
|
2776 |
|
2777 // I don't think there's any other way to revert back to unknown |
|
2778 // constant types, so here it is. |
|
2779 |
|
2780 if (old_type_tag == unknown_constant && error_state) |
|
2781 { |
|
2782 if (type_tag == matrix_constant) |
|
2783 delete matrix; |
|
2784 else if (type_tag == complex_matrix_constant) |
|
2785 delete complex_matrix; |
|
2786 |
|
2787 type_tag = unknown_constant; |
|
2788 } |
|
2789 } |
|
2790 else if (nargs > 3 || nargs < 2) |
|
2791 ::error ("invalid index expression for scalar type"); |
|
2792 else |
|
2793 ::error ("index invalid or out of range for scalar type"); |
|
2794 } |
|
2795 |
|
2796 /* |
|
2797 * Assignments to matrices (and vectors). |
|
2798 * |
|
2799 * For compatibility with Matlab, we allow assignment of an empty |
|
2800 * matrix to an expression with empty indices to do nothing. |
|
2801 */ |
|
2802 void |
|
2803 tree_constant_rep::do_matrix_assignment (tree_constant& rhs, |
500
|
2804 const Octave_object& args, int nargs) |
492
|
2805 { |
|
2806 assert (type_tag == unknown_constant |
|
2807 || type_tag == matrix_constant |
|
2808 || type_tag == complex_matrix_constant); |
|
2809 |
|
2810 if (type_tag == matrix_constant && rhs.is_complex_type ()) |
|
2811 { |
|
2812 Matrix *old_matrix = matrix; |
|
2813 complex_matrix = new ComplexMatrix (*matrix); |
|
2814 type_tag = complex_matrix_constant; |
|
2815 delete old_matrix; |
|
2816 } |
|
2817 else if (type_tag == unknown_constant) |
|
2818 { |
|
2819 if (rhs.is_complex_type ()) |
|
2820 { |
|
2821 complex_matrix = new ComplexMatrix (); |
|
2822 type_tag = complex_matrix_constant; |
|
2823 } |
|
2824 else |
|
2825 { |
|
2826 matrix = new Matrix (); |
|
2827 type_tag = matrix_constant; |
|
2828 } |
|
2829 } |
|
2830 |
|
2831 // The do_matrix_assignment functions can't handle empty matrices, so |
|
2832 // don't let any pass through here. |
|
2833 switch (nargs) |
|
2834 { |
|
2835 case 2: |
500
|
2836 if (args.length () <= 0) |
492
|
2837 ::error ("matrix index is null"); |
500
|
2838 else if (args(1).is_undefined ()) |
492
|
2839 ::error ("matrix index is undefined"); |
|
2840 else |
500
|
2841 do_matrix_assignment (rhs, args(1)); |
492
|
2842 break; |
|
2843 case 3: |
500
|
2844 if (args.length () <= 0) |
492
|
2845 ::error ("matrix indices are null"); |
500
|
2846 else if (args(1).is_undefined ()) |
492
|
2847 ::error ("first matrix index is undefined"); |
500
|
2848 else if (args(2).is_undefined ()) |
492
|
2849 ::error ("second matrix index is undefined"); |
500
|
2850 else if (args(1).is_empty () || args(2).is_empty ()) |
492
|
2851 { |
|
2852 if (! rhs.is_empty ()) |
|
2853 { |
|
2854 ::error ("in assignment expression, a matrix index is empty"); |
|
2855 ::error ("but hte right hand side is not an empty matrix"); |
|
2856 } |
|
2857 // XXX FIXME XXX -- to really be correct here, we should probably |
|
2858 // check to see if the assignment conforms, but that seems like more |
|
2859 // work than it's worth right now... |
|
2860 } |
|
2861 else |
500
|
2862 do_matrix_assignment (rhs, args(1), args(2)); |
492
|
2863 break; |
|
2864 default: |
|
2865 ::error ("too many indices for matrix expression"); |
|
2866 break; |
|
2867 } |
|
2868 } |
|
2869 |
|
2870 /* |
|
2871 * Matrix assignments indexed by a single value. |
|
2872 */ |
|
2873 void |
|
2874 tree_constant_rep::do_matrix_assignment (tree_constant& rhs, |
|
2875 tree_constant& i_arg) |
|
2876 { |
|
2877 int nr = rows (); |
|
2878 int nc = columns (); |
|
2879 |
|
2880 if (user_pref.do_fortran_indexing || nr <= 1 || nc <= 1) |
|
2881 { |
|
2882 if (i_arg.is_empty ()) |
|
2883 { |
|
2884 if (! rhs.is_empty ()) |
|
2885 { |
|
2886 ::error ("in assignment expression, matrix index is empty but"); |
|
2887 ::error ("right hand side is not an empty matrix"); |
|
2888 } |
|
2889 // XXX FIXME XXX -- to really be correct here, we should probably |
|
2890 // check to see if the assignment conforms, but that seems like more |
|
2891 // work than it's worth right now... |
|
2892 |
|
2893 // The assignment functions can't handle empty matrices, so don't let |
|
2894 // any pass through here. |
|
2895 return; |
|
2896 } |
|
2897 |
|
2898 // We can't handle the case of assigning to a vector first, since even |
|
2899 // then, the two operations are not equivalent. For example, the |
|
2900 // expression V(:) = M is handled differently depending on whether the |
|
2901 // user specified do_fortran_indexing = "true". |
|
2902 |
|
2903 if (user_pref.do_fortran_indexing) |
|
2904 fortran_style_matrix_assignment (rhs, i_arg); |
|
2905 else if (nr <= 1 || nc <= 1) |
|
2906 vector_assignment (rhs, i_arg); |
|
2907 else |
|
2908 panic_impossible (); |
|
2909 } |
|
2910 else |
|
2911 ::error ("single index only valid for row or column vector"); |
|
2912 } |
|
2913 |
|
2914 /* |
|
2915 * Fortran-style assignments. Matrices are assumed to be stored in |
|
2916 * column-major order and it is ok to use a single index for |
|
2917 * multi-dimensional matrices. |
|
2918 */ |
|
2919 void |
|
2920 tree_constant_rep::fortran_style_matrix_assignment (tree_constant& rhs, |
|
2921 tree_constant& i_arg) |
|
2922 { |
|
2923 tree_constant tmp_i = i_arg.make_numeric_or_magic (); |
|
2924 |
|
2925 tree_constant_rep::constant_type itype = tmp_i.const_type (); |
|
2926 |
|
2927 int nr = rows (); |
|
2928 int nc = columns (); |
|
2929 |
|
2930 int rhs_nr = rhs.rows (); |
|
2931 int rhs_nc = rhs.columns (); |
|
2932 |
|
2933 switch (itype) |
|
2934 { |
|
2935 case complex_scalar_constant: |
|
2936 case scalar_constant: |
|
2937 { |
|
2938 int i = NINT (tmp_i.double_value ()); |
|
2939 int idx = i - 1; |
|
2940 |
|
2941 if (rhs_nr == 0 && rhs_nc == 0) |
|
2942 { |
|
2943 if (idx < nr * nc) |
|
2944 { |
|
2945 convert_to_row_or_column_vector (); |
|
2946 |
|
2947 nr = rows (); |
|
2948 nc = columns (); |
|
2949 |
|
2950 if (nr == 1) |
|
2951 delete_column (idx); |
|
2952 else if (nc == 1) |
|
2953 delete_row (idx); |
|
2954 else |
|
2955 panic_impossible (); |
|
2956 } |
|
2957 return; |
|
2958 } |
|
2959 |
|
2960 if (index_check (idx, "") < 0) |
|
2961 return; |
|
2962 |
|
2963 if (nr <= 1 || nc <= 1) |
|
2964 { |
|
2965 maybe_resize (idx); |
|
2966 if (error_state) |
|
2967 return; |
|
2968 } |
|
2969 else if (range_max_check (idx, nr * nc) < 0) |
|
2970 return; |
|
2971 |
|
2972 nr = rows (); |
|
2973 nc = columns (); |
|
2974 |
|
2975 if (! indexed_assign_conforms (1, 1, rhs_nr, rhs_nc)) |
|
2976 { |
|
2977 ::error ("for A(int) = X: X must be a scalar"); |
|
2978 return; |
|
2979 } |
|
2980 int ii = fortran_row (i, nr) - 1; |
|
2981 int jj = fortran_column (i, nr) - 1; |
|
2982 do_matrix_assignment (rhs, ii, jj); |
|
2983 } |
|
2984 break; |
|
2985 case complex_matrix_constant: |
|
2986 case matrix_constant: |
|
2987 { |
|
2988 Matrix mi = tmp_i.matrix_value (); |
|
2989 int len = nr * nc; |
|
2990 idx_vector ii (mi, 1, "", len); // Always do fortran indexing here... |
|
2991 if (! ii) |
|
2992 return; |
|
2993 |
|
2994 if (rhs_nr == 0 && rhs_nc == 0) |
|
2995 { |
|
2996 ii.sort_uniq (); |
|
2997 int num_to_delete = 0; |
|
2998 for (int i = 0; i < ii.length (); i++) |
|
2999 { |
|
3000 if (ii.elem (i) < len) |
|
3001 num_to_delete++; |
|
3002 else |
|
3003 break; |
|
3004 } |
|
3005 |
|
3006 if (num_to_delete > 0) |
|
3007 { |
|
3008 if (num_to_delete != ii.length ()) |
|
3009 ii.shorten (num_to_delete); |
|
3010 |
|
3011 convert_to_row_or_column_vector (); |
|
3012 |
|
3013 nr = rows (); |
|
3014 nc = columns (); |
|
3015 |
|
3016 if (nr == 1) |
|
3017 delete_columns (ii); |
|
3018 else if (nc == 1) |
|
3019 delete_rows (ii); |
|
3020 else |
|
3021 panic_impossible (); |
|
3022 } |
|
3023 return; |
|
3024 } |
|
3025 |
|
3026 if (nr <= 1 || nc <= 1) |
|
3027 { |
|
3028 maybe_resize (ii.max ()); |
|
3029 if (error_state) |
|
3030 return; |
|
3031 } |
|
3032 else if (range_max_check (ii.max (), len) < 0) |
|
3033 return; |
|
3034 |
|
3035 int ilen = ii.capacity (); |
|
3036 |
|
3037 if (ilen != rhs_nr * rhs_nc) |
|
3038 { |
|
3039 ::error ("A(matrix) = X: X and matrix must have the same number"); |
|
3040 ::error ("of elements"); |
|
3041 } |
|
3042 else if (ilen == 1 && rhs.is_scalar_type ()) |
|
3043 { |
|
3044 int nr = rows (); |
|
3045 int idx = ii.elem (0); |
|
3046 int ii = fortran_row (idx + 1, nr) - 1; |
|
3047 int jj = fortran_column (idx + 1, nr) - 1; |
|
3048 |
|
3049 if (rhs.const_type () == scalar_constant) |
|
3050 matrix->elem (ii, jj) = rhs.double_value (); |
|
3051 else if (rhs.const_type () == complex_scalar_constant) |
|
3052 complex_matrix->elem (ii, jj) = rhs.complex_value (); |
|
3053 else |
|
3054 panic_impossible (); |
|
3055 } |
|
3056 else |
|
3057 fortran_style_matrix_assignment (rhs, ii); |
|
3058 } |
|
3059 break; |
|
3060 case string_constant: |
|
3061 gripe_string_invalid (); |
|
3062 break; |
|
3063 case range_constant: |
|
3064 gripe_range_invalid (); |
|
3065 break; |
|
3066 case magic_colon: |
|
3067 // a(:) = [] is equivalent to a(:,:) = []. |
|
3068 if (rhs_nr == 0 && rhs_nc == 0) |
|
3069 do_matrix_assignment (rhs, magic_colon, magic_colon); |
|
3070 else |
|
3071 fortran_style_matrix_assignment (rhs, magic_colon); |
|
3072 break; |
|
3073 default: |
|
3074 panic_impossible (); |
|
3075 break; |
|
3076 } |
|
3077 } |
|
3078 |
|
3079 /* |
|
3080 * Fortran-style assignment for vector index. |
|
3081 */ |
|
3082 void |
|
3083 tree_constant_rep::fortran_style_matrix_assignment (tree_constant& rhs, |
|
3084 idx_vector& i) |
|
3085 { |
|
3086 assert (rhs.is_matrix_type ()); |
|
3087 |
|
3088 int ilen = i.capacity (); |
|
3089 |
|
3090 REP_RHS_MATRIX (rhs, rhs_m, rhs_cm, rhs_nr, rhs_nc); |
|
3091 |
|
3092 int len = rhs_nr * rhs_nc; |
|
3093 |
|
3094 if (len == ilen) |
|
3095 { |
|
3096 int nr = rows (); |
|
3097 if (rhs.const_type () == matrix_constant) |
|
3098 { |
|
3099 double *cop_out = rhs_m.fortran_vec (); |
|
3100 for (int k = 0; k < len; k++) |
|
3101 { |
|
3102 int ii = fortran_row (i.elem (k) + 1, nr) - 1; |
|
3103 int jj = fortran_column (i.elem (k) + 1, nr) - 1; |
|
3104 |
|
3105 matrix->elem (ii, jj) = *cop_out++; |
|
3106 } |
|
3107 } |
|
3108 else |
|
3109 { |
|
3110 Complex *cop_out = rhs_cm.fortran_vec (); |
|
3111 for (int k = 0; k < len; k++) |
|
3112 { |
|
3113 int ii = fortran_row (i.elem (k) + 1, nr) - 1; |
|
3114 int jj = fortran_column (i.elem (k) + 1, nr) - 1; |
|
3115 |
|
3116 complex_matrix->elem (ii, jj) = *cop_out++; |
|
3117 } |
|
3118 } |
|
3119 } |
|
3120 else |
|
3121 ::error ("number of rows and columns must match for indexed assignment"); |
|
3122 } |
|
3123 |
|
3124 /* |
|
3125 * Fortran-style assignment for colon index. |
|
3126 */ |
|
3127 void |
|
3128 tree_constant_rep::fortran_style_matrix_assignment |
|
3129 (tree_constant& rhs, tree_constant_rep::constant_type mci) |
|
3130 { |
|
3131 assert (rhs.is_matrix_type () && mci == tree_constant_rep::magic_colon); |
|
3132 |
|
3133 int nr = rows (); |
|
3134 int nc = columns (); |
|
3135 |
|
3136 REP_RHS_MATRIX (rhs, rhs_m, rhs_cm, rhs_nr, rhs_nc); |
|
3137 |
|
3138 int rhs_size = rhs_nr * rhs_nc; |
|
3139 if (rhs_size == 0) |
|
3140 { |
|
3141 if (rhs.const_type () == matrix_constant) |
|
3142 { |
|
3143 delete matrix; |
|
3144 matrix = new Matrix (0, 0); |
|
3145 return; |
|
3146 } |
|
3147 else |
|
3148 panic_impossible (); |
|
3149 } |
|
3150 else if (nr*nc != rhs_size) |
|
3151 { |
|
3152 ::error ("A(:) = X: X and A must have the same number of elements"); |
|
3153 return; |
|
3154 } |
|
3155 |
|
3156 if (rhs.const_type () == matrix_constant) |
|
3157 { |
|
3158 double *cop_out = rhs_m.fortran_vec (); |
|
3159 for (int j = 0; j < nc; j++) |
|
3160 for (int i = 0; i < nr; i++) |
|
3161 matrix->elem (i, j) = *cop_out++; |
|
3162 } |
|
3163 else |
|
3164 { |
|
3165 Complex *cop_out = rhs_cm.fortran_vec (); |
|
3166 for (int j = 0; j < nc; j++) |
|
3167 for (int i = 0; i < nr; i++) |
|
3168 complex_matrix->elem (i, j) = *cop_out++; |
|
3169 } |
|
3170 } |
|
3171 |
|
3172 /* |
|
3173 * Assignments to vectors. Hand off to other functions once we know |
|
3174 * what kind of index we have. For a colon, it is the same as |
|
3175 * assignment to a matrix indexed by two colons. |
|
3176 */ |
|
3177 void |
|
3178 tree_constant_rep::vector_assignment (tree_constant& rhs, tree_constant& i_arg) |
|
3179 { |
|
3180 int nr = rows (); |
|
3181 int nc = columns (); |
|
3182 |
|
3183 assert ((nr == 1 || nc == 1 || (nr == 0 && nc == 0)) |
|
3184 && ! user_pref.do_fortran_indexing); |
|
3185 |
|
3186 tree_constant tmp_i = i_arg.make_numeric_or_range_or_magic (); |
|
3187 |
|
3188 tree_constant_rep::constant_type itype = tmp_i.const_type (); |
|
3189 |
|
3190 switch (itype) |
|
3191 { |
|
3192 case complex_scalar_constant: |
|
3193 case scalar_constant: |
|
3194 { |
|
3195 int i = tree_to_mat_idx (tmp_i.double_value ()); |
|
3196 if (index_check (i, "") < 0) |
|
3197 return; |
|
3198 do_vector_assign (rhs, i); |
|
3199 } |
|
3200 break; |
|
3201 case complex_matrix_constant: |
|
3202 case matrix_constant: |
|
3203 { |
|
3204 Matrix mi = tmp_i.matrix_value (); |
|
3205 int len = nr * nc; |
|
3206 idx_vector iv (mi, user_pref.do_fortran_indexing, "", len); |
|
3207 if (! iv) |
|
3208 return; |
|
3209 |
|
3210 do_vector_assign (rhs, iv); |
|
3211 } |
|
3212 break; |
|
3213 case string_constant: |
|
3214 gripe_string_invalid (); |
|
3215 break; |
|
3216 case range_constant: |
|
3217 { |
|
3218 Range ri = tmp_i.range_value (); |
|
3219 int len = nr * nc; |
|
3220 if (len == 2 && is_zero_one (ri)) |
|
3221 { |
|
3222 do_vector_assign (rhs, 1); |
|
3223 } |
|
3224 else if (len == 2 && is_one_zero (ri)) |
|
3225 { |
|
3226 do_vector_assign (rhs, 0); |
|
3227 } |
|
3228 else |
|
3229 { |
|
3230 if (index_check (ri, "") < 0) |
|
3231 return; |
|
3232 do_vector_assign (rhs, ri); |
|
3233 } |
|
3234 } |
|
3235 break; |
|
3236 case magic_colon: |
|
3237 { |
|
3238 int rhs_nr = rhs.rows (); |
|
3239 int rhs_nc = rhs.columns (); |
|
3240 |
|
3241 if (! indexed_assign_conforms (nr, nc, rhs_nr, rhs_nc)) |
|
3242 { |
|
3243 ::error ("A(:) = X: X and A must have the same dimensions"); |
|
3244 return; |
|
3245 } |
|
3246 do_matrix_assignment (rhs, magic_colon, magic_colon); |
|
3247 } |
|
3248 break; |
|
3249 default: |
|
3250 panic_impossible (); |
|
3251 break; |
|
3252 } |
|
3253 } |
|
3254 |
|
3255 /* |
|
3256 * Check whether an indexed assignment to a vector is valid. |
|
3257 */ |
|
3258 void |
|
3259 tree_constant_rep::check_vector_assign (int rhs_nr, int rhs_nc, |
|
3260 int ilen, const char *rm) |
|
3261 { |
|
3262 int nr = rows (); |
|
3263 int nc = columns (); |
|
3264 |
|
3265 if ((nr == 1 && nc == 1) || nr == 0 || nc == 0) // No orientation. |
|
3266 { |
|
3267 if (! (ilen == rhs_nr || ilen == rhs_nc)) |
|
3268 { |
|
3269 ::error ("A(%s) = X: X and %s must have the same number of elements", |
|
3270 rm, rm); |
|
3271 } |
|
3272 } |
|
3273 else if (nr == 1) // Preserve current row orientation. |
|
3274 { |
|
3275 if (! (rhs_nr == 1 && rhs_nc == ilen)) |
|
3276 { |
|
3277 ::error ("A(%s) = X: where A is a row vector, X must also be a", rm); |
|
3278 ::error ("row vector with the same number of elements as %s", rm); |
|
3279 } |
|
3280 } |
|
3281 else if (nc == 1) // Preserve current column orientation. |
|
3282 { |
|
3283 if (! (rhs_nc == 1 && rhs_nr == ilen)) |
|
3284 { |
|
3285 ::error ("A(%s) = X: where A is a column vector, X must also be", rm); |
|
3286 ::error ("a column vector with the same number of elements as %s", rm); |
|
3287 } |
|
3288 } |
|
3289 else |
|
3290 panic_impossible (); |
|
3291 } |
|
3292 |
|
3293 /* |
|
3294 * Assignment to a vector with an integer index. |
|
3295 */ |
|
3296 void |
|
3297 tree_constant_rep::do_vector_assign (tree_constant& rhs, int i) |
|
3298 { |
|
3299 int rhs_nr = rhs.rows (); |
|
3300 int rhs_nc = rhs.columns (); |
|
3301 |
|
3302 if (indexed_assign_conforms (1, 1, rhs_nr, rhs_nc)) |
|
3303 { |
|
3304 maybe_resize (i); |
|
3305 if (error_state) |
|
3306 return; |
|
3307 |
|
3308 int nr = rows (); |
|
3309 int nc = columns (); |
|
3310 |
|
3311 if (nr == 1) |
|
3312 { |
|
3313 REP_ELEM_ASSIGN (0, i, rhs.double_value (), rhs.complex_value (), |
|
3314 rhs.is_real_type ()); |
|
3315 } |
|
3316 else if (nc == 1) |
|
3317 { |
|
3318 REP_ELEM_ASSIGN (i, 0, rhs.double_value (), rhs.complex_value (), |
|
3319 rhs.is_real_type ()); |
|
3320 } |
|
3321 else |
|
3322 panic_impossible (); |
|
3323 } |
|
3324 else if (rhs_nr == 0 && rhs_nc == 0) |
|
3325 { |
|
3326 int nr = rows (); |
|
3327 int nc = columns (); |
|
3328 |
|
3329 int len = MAX (nr, nc); |
|
3330 |
|
3331 if (i < 0 || i >= len) |
|
3332 { |
|
3333 ::error ("A(int) = []: index out of range"); |
|
3334 return; |
|
3335 } |
|
3336 |
|
3337 if (nr == 1) |
|
3338 delete_column (i); |
|
3339 else if (nc == 1) |
|
3340 delete_row (i); |
|
3341 else |
|
3342 panic_impossible (); |
|
3343 } |
|
3344 else |
|
3345 { |
|
3346 ::error ("for A(int) = X: X must be a scalar"); |
|
3347 return; |
|
3348 } |
|
3349 } |
|
3350 |
|
3351 /* |
|
3352 * Assignment to a vector with a vector index. |
|
3353 */ |
|
3354 void |
|
3355 tree_constant_rep::do_vector_assign (tree_constant& rhs, idx_vector& iv) |
|
3356 { |
|
3357 if (rhs.is_zero_by_zero ()) |
|
3358 { |
|
3359 int nr = rows (); |
|
3360 int nc = columns (); |
|
3361 |
|
3362 int len = MAX (nr, nc); |
|
3363 |
|
3364 if (iv.max () >= len) |
|
3365 { |
|
3366 ::error ("A(matrix) = []: index out of range"); |
|
3367 return; |
|
3368 } |
|
3369 |
|
3370 if (nr == 1) |
|
3371 delete_columns (iv); |
|
3372 else if (nc == 1) |
|
3373 delete_rows (iv); |
|
3374 else |
|
3375 panic_impossible (); |
|
3376 } |
|
3377 else if (rhs.is_scalar_type ()) |
|
3378 { |
|
3379 int nr = rows (); |
|
3380 int nc = columns (); |
|
3381 |
|
3382 if (iv.capacity () == 1) |
|
3383 { |
|
3384 int idx = iv.elem (0); |
|
3385 |
|
3386 if (nr == 1) |
|
3387 { |
|
3388 REP_ELEM_ASSIGN (0, idx, rhs.double_value (), |
|
3389 rhs.complex_value (), rhs.is_real_type ()); |
|
3390 } |
|
3391 else if (nc == 1) |
|
3392 { |
|
3393 REP_ELEM_ASSIGN (idx, 0, rhs.double_value (), |
|
3394 rhs.complex_value (), rhs.is_real_type ()); |
|
3395 } |
|
3396 else |
|
3397 panic_impossible (); |
|
3398 } |
|
3399 else |
|
3400 { |
|
3401 if (nr == 1) |
|
3402 { |
|
3403 ::error ("A(matrix) = X: where A is a row vector, X must also be a"); |
|
3404 ::error ("row vector with the same number of elements as matrix"); |
|
3405 } |
|
3406 else if (nc == 1) |
|
3407 { |
|
3408 ::error ("A(matrix) = X: where A is a column vector, X must also be a"); |
|
3409 ::error ("column vector with the same number of elements as matrix"); |
|
3410 } |
|
3411 else |
|
3412 panic_impossible (); |
|
3413 } |
|
3414 } |
|
3415 else if (rhs.is_matrix_type ()) |
|
3416 { |
|
3417 REP_RHS_MATRIX (rhs, rhs_m, rhs_cm, rhs_nr, rhs_nc); |
|
3418 |
|
3419 int ilen = iv.capacity (); |
|
3420 check_vector_assign (rhs_nr, rhs_nc, ilen, "matrix"); |
|
3421 if (error_state) |
|
3422 return; |
|
3423 |
|
3424 force_orient f_orient = no_orient; |
|
3425 if (rhs_nr == 1 && rhs_nc != 1) |
|
3426 f_orient = row_orient; |
|
3427 else if (rhs_nc == 1 && rhs_nr != 1) |
|
3428 f_orient = column_orient; |
|
3429 |
|
3430 maybe_resize (iv.max (), f_orient); |
|
3431 if (error_state) |
|
3432 return; |
|
3433 |
|
3434 int nr = rows (); |
|
3435 int nc = columns (); |
|
3436 |
|
3437 if (nr == 1) |
|
3438 { |
|
3439 for (int i = 0; i < iv.capacity (); i++) |
|
3440 REP_ELEM_ASSIGN (0, iv.elem (i), rhs_m.elem (0, i), |
|
3441 rhs_cm.elem (0, i), rhs.is_real_type ()); |
|
3442 } |
|
3443 else if (nc == 1) |
|
3444 { |
|
3445 for (int i = 0; i < iv.capacity (); i++) |
|
3446 REP_ELEM_ASSIGN (iv.elem (i), 0, rhs_m.elem (i, 0), |
|
3447 rhs_cm.elem (i, 0), rhs.is_real_type ()); |
|
3448 } |
|
3449 else |
|
3450 panic_impossible (); |
|
3451 } |
|
3452 else |
|
3453 panic_impossible (); |
|
3454 } |
|
3455 |
|
3456 /* |
|
3457 * Assignment to a vector with a range index. |
|
3458 */ |
|
3459 void |
|
3460 tree_constant_rep::do_vector_assign (tree_constant& rhs, Range& ri) |
|
3461 { |
|
3462 if (rhs.is_zero_by_zero ()) |
|
3463 { |
|
3464 int nr = rows (); |
|
3465 int nc = columns (); |
|
3466 |
|
3467 int len = MAX (nr, nc); |
|
3468 |
|
3469 int b = tree_to_mat_idx (ri.min ()); |
|
3470 int l = tree_to_mat_idx (ri.max ()); |
|
3471 if (b < 0 || l >= len) |
|
3472 { |
|
3473 ::error ("A(range) = []: index out of range"); |
|
3474 return; |
|
3475 } |
|
3476 |
|
3477 if (nr == 1) |
|
3478 delete_columns (ri); |
|
3479 else if (nc == 1) |
|
3480 delete_rows (ri); |
|
3481 else |
|
3482 panic_impossible (); |
|
3483 } |
|
3484 else if (rhs.is_scalar_type ()) |
|
3485 { |
|
3486 int nr = rows (); |
|
3487 int nc = columns (); |
|
3488 |
|
3489 if (nr == 1) |
|
3490 { |
|
3491 ::error ("A(range) = X: where A is a row vector, X must also be a"); |
|
3492 ::error ("row vector with the same number of elements as range"); |
|
3493 } |
|
3494 else if (nc == 1) |
|
3495 { |
|
3496 ::error ("A(range) = X: where A is a column vector, X must also be a"); |
|
3497 ::error ("column vector with the same number of elements as range"); |
|
3498 } |
|
3499 else |
|
3500 panic_impossible (); |
|
3501 } |
|
3502 else if (rhs.is_matrix_type ()) |
|
3503 { |
|
3504 REP_RHS_MATRIX (rhs, rhs_m, rhs_cm, rhs_nr, rhs_nc); |
|
3505 |
|
3506 int ilen = ri.nelem (); |
|
3507 check_vector_assign (rhs_nr, rhs_nc, ilen, "range"); |
|
3508 if (error_state) |
|
3509 return; |
|
3510 |
|
3511 force_orient f_orient = no_orient; |
|
3512 if (rhs_nr == 1 && rhs_nc != 1) |
|
3513 f_orient = row_orient; |
|
3514 else if (rhs_nc == 1 && rhs_nr != 1) |
|
3515 f_orient = column_orient; |
|
3516 |
|
3517 maybe_resize (tree_to_mat_idx (ri.max ()), f_orient); |
|
3518 if (error_state) |
|
3519 return; |
|
3520 |
|
3521 int nr = rows (); |
|
3522 int nc = columns (); |
|
3523 |
|
3524 double b = ri.base (); |
|
3525 double increment = ri.inc (); |
|
3526 |
|
3527 if (nr == 1) |
|
3528 { |
|
3529 for (int i = 0; i < ri.nelem (); i++) |
|
3530 { |
|
3531 double tmp = b + i * increment; |
|
3532 int col = tree_to_mat_idx (tmp); |
|
3533 REP_ELEM_ASSIGN (0, col, rhs_m.elem (0, i), rhs_cm.elem (0, i), |
|
3534 rhs.is_real_type ()); |
|
3535 } |
|
3536 } |
|
3537 else if (nc == 1) |
|
3538 { |
|
3539 for (int i = 0; i < ri.nelem (); i++) |
|
3540 { |
|
3541 double tmp = b + i * increment; |
|
3542 int row = tree_to_mat_idx (tmp); |
|
3543 REP_ELEM_ASSIGN (row, 0, rhs_m.elem (i, 0), rhs_cm.elem (i, 0), |
|
3544 rhs.is_real_type ()); |
|
3545 } |
|
3546 } |
|
3547 else |
|
3548 panic_impossible (); |
|
3549 } |
|
3550 else |
|
3551 panic_impossible (); |
|
3552 } |
|
3553 |
|
3554 /* |
|
3555 * Matrix assignment indexed by two values. This function determines |
|
3556 * the type of the first arugment, checks as much as possible, and |
|
3557 * then calls one of a set of functions to handle the specific cases: |
|
3558 * |
|
3559 * M (integer, arg2) = RHS (MA1) |
|
3560 * M (vector, arg2) = RHS (MA2) |
|
3561 * M (range, arg2) = RHS (MA3) |
|
3562 * M (colon, arg2) = RHS (MA4) |
|
3563 * |
|
3564 * Each of those functions determines the type of the second argument |
|
3565 * and calls another function to handle the real work of doing the |
|
3566 * assignment. |
|
3567 */ |
|
3568 void |
|
3569 tree_constant_rep::do_matrix_assignment (tree_constant& rhs, |
|
3570 tree_constant& i_arg, |
|
3571 tree_constant& j_arg) |
|
3572 { |
|
3573 tree_constant tmp_i = i_arg.make_numeric_or_range_or_magic (); |
|
3574 |
|
3575 tree_constant_rep::constant_type itype = tmp_i.const_type (); |
|
3576 |
|
3577 switch (itype) |
|
3578 { |
|
3579 case complex_scalar_constant: |
|
3580 case scalar_constant: |
|
3581 { |
|
3582 int i = tree_to_mat_idx (tmp_i.double_value ()); |
|
3583 if (index_check (i, "row") < 0) |
|
3584 return; |
|
3585 do_matrix_assignment (rhs, i, j_arg); |
|
3586 } |
|
3587 break; |
|
3588 case complex_matrix_constant: |
|
3589 case matrix_constant: |
|
3590 { |
|
3591 Matrix mi = tmp_i.matrix_value (); |
|
3592 idx_vector iv (mi, user_pref.do_fortran_indexing, "row", rows ()); |
|
3593 if (! iv) |
|
3594 return; |
|
3595 |
|
3596 do_matrix_assignment (rhs, iv, j_arg); |
|
3597 } |
|
3598 break; |
|
3599 case string_constant: |
|
3600 gripe_string_invalid (); |
|
3601 break; |
|
3602 case range_constant: |
|
3603 { |
|
3604 Range ri = tmp_i.range_value (); |
|
3605 int nr = rows (); |
|
3606 if (nr == 2 && is_zero_one (ri)) |
|
3607 { |
|
3608 do_matrix_assignment (rhs, 1, j_arg); |
|
3609 } |
|
3610 else if (nr == 2 && is_one_zero (ri)) |
|
3611 { |
|
3612 do_matrix_assignment (rhs, 0, j_arg); |
|
3613 } |
|
3614 else |
|
3615 { |
|
3616 if (index_check (ri, "row") < 0) |
|
3617 return; |
|
3618 do_matrix_assignment (rhs, ri, j_arg); |
|
3619 } |
|
3620 } |
|
3621 break; |
|
3622 case magic_colon: |
|
3623 do_matrix_assignment (rhs, magic_colon, j_arg); |
|
3624 break; |
|
3625 default: |
|
3626 panic_impossible (); |
|
3627 break; |
|
3628 } |
|
3629 } |
|
3630 |
|
3631 /* MA1 */ |
|
3632 void |
|
3633 tree_constant_rep::do_matrix_assignment (tree_constant& rhs, int i, |
|
3634 tree_constant& j_arg) |
|
3635 { |
|
3636 tree_constant tmp_j = j_arg.make_numeric_or_range_or_magic (); |
|
3637 |
|
3638 tree_constant_rep::constant_type jtype = tmp_j.const_type (); |
|
3639 |
|
3640 int rhs_nr = rhs.rows (); |
|
3641 int rhs_nc = rhs.columns (); |
|
3642 |
|
3643 switch (jtype) |
|
3644 { |
|
3645 case complex_scalar_constant: |
|
3646 case scalar_constant: |
|
3647 { |
|
3648 int j = tree_to_mat_idx (tmp_j.double_value ()); |
|
3649 if (index_check (j, "column") < 0) |
|
3650 return; |
|
3651 if (! indexed_assign_conforms (1, 1, rhs_nr, rhs_nc)) |
|
3652 { |
|
3653 ::error ("A(int,int) = X, X must be a scalar"); |
|
3654 return; |
|
3655 } |
|
3656 maybe_resize (i, j); |
|
3657 if (error_state) |
|
3658 return; |
|
3659 |
|
3660 do_matrix_assignment (rhs, i, j); |
|
3661 } |
|
3662 break; |
|
3663 case complex_matrix_constant: |
|
3664 case matrix_constant: |
|
3665 { |
|
3666 Matrix mj = tmp_j.matrix_value (); |
|
3667 idx_vector jv (mj, user_pref.do_fortran_indexing, "column", |
|
3668 columns ()); |
|
3669 if (! jv) |
|
3670 return; |
|
3671 |
|
3672 if (! indexed_assign_conforms (1, jv.capacity (), rhs_nr, rhs_nc)) |
|
3673 { |
|
3674 ::error ("A(int,matrix) = X: X must be a row vector with the same"); |
|
3675 ::error ("number of elements as matrix"); |
|
3676 return; |
|
3677 } |
|
3678 maybe_resize (i, jv.max ()); |
|
3679 if (error_state) |
|
3680 return; |
|
3681 |
|
3682 do_matrix_assignment (rhs, i, jv); |
|
3683 } |
|
3684 break; |
|
3685 case string_constant: |
|
3686 gripe_string_invalid (); |
|
3687 break; |
|
3688 case range_constant: |
|
3689 { |
|
3690 Range rj = tmp_j.range_value (); |
|
3691 if (! indexed_assign_conforms (1, rj.nelem (), rhs_nr, rhs_nc)) |
|
3692 { |
|
3693 ::error ("A(int,range) = X: X must be a row vector with the same"); |
|
3694 ::error ("number of elements as range"); |
|
3695 return; |
|
3696 } |
|
3697 |
|
3698 int nc = columns (); |
|
3699 if (nc == 2 && is_zero_one (rj) && rhs_nc == 1) |
|
3700 { |
|
3701 do_matrix_assignment (rhs, i, 1); |
|
3702 } |
|
3703 else if (nc == 2 && is_one_zero (rj) && rhs_nc == 1) |
|
3704 { |
|
3705 do_matrix_assignment (rhs, i, 0); |
|
3706 } |
|
3707 else |
|
3708 { |
|
3709 if (index_check (rj, "column") < 0) |
|
3710 return; |
|
3711 maybe_resize (i, tree_to_mat_idx (rj.max ())); |
|
3712 if (error_state) |
|
3713 return; |
|
3714 |
|
3715 do_matrix_assignment (rhs, i, rj); |
|
3716 } |
|
3717 } |
|
3718 break; |
|
3719 case magic_colon: |
|
3720 { |
|
3721 int nc = columns (); |
|
3722 int nr = rows (); |
|
3723 if (nc == 0 && nr == 0 && rhs_nr == 1) |
|
3724 { |
|
3725 if (rhs.is_complex_type ()) |
|
3726 { |
|
3727 complex_matrix = new ComplexMatrix (); |
|
3728 type_tag = complex_matrix_constant; |
|
3729 } |
|
3730 else |
|
3731 { |
|
3732 matrix = new Matrix (); |
|
3733 type_tag = matrix_constant; |
|
3734 } |
|
3735 maybe_resize (i, rhs_nc-1); |
|
3736 if (error_state) |
|
3737 return; |
|
3738 } |
|
3739 else if (indexed_assign_conforms (1, nc, rhs_nr, rhs_nc)) |
|
3740 { |
|
3741 maybe_resize (i, nc-1); |
|
3742 if (error_state) |
|
3743 return; |
|
3744 } |
|
3745 else if (rhs_nr == 0 && rhs_nc == 0) |
|
3746 { |
|
3747 if (i < 0 || i >= nr) |
|
3748 { |
|
3749 ::error ("A(int,:) = []: row index out of range"); |
|
3750 return; |
|
3751 } |
|
3752 } |
|
3753 else |
|
3754 { |
|
3755 ::error ("A(int,:) = X: X must be a row vector with the same"); |
|
3756 ::error ("number of columns as A"); |
|
3757 return; |
|
3758 } |
|
3759 |
|
3760 do_matrix_assignment (rhs, i, magic_colon); |
|
3761 } |
|
3762 break; |
|
3763 default: |
|
3764 panic_impossible (); |
|
3765 break; |
|
3766 } |
|
3767 } |
|
3768 |
|
3769 /* MA2 */ |
|
3770 void |
|
3771 tree_constant_rep::do_matrix_assignment (tree_constant& rhs, idx_vector& iv, |
|
3772 tree_constant& j_arg) |
|
3773 { |
|
3774 tree_constant tmp_j = j_arg.make_numeric_or_range_or_magic (); |
|
3775 |
|
3776 tree_constant_rep::constant_type jtype = tmp_j.const_type (); |
|
3777 |
|
3778 int rhs_nr = rhs.rows (); |
|
3779 int rhs_nc = rhs.columns (); |
|
3780 |
|
3781 switch (jtype) |
|
3782 { |
|
3783 case complex_scalar_constant: |
|
3784 case scalar_constant: |
|
3785 { |
|
3786 int j = tree_to_mat_idx (tmp_j.double_value ()); |
|
3787 if (index_check (j, "column") < 0) |
|
3788 return; |
|
3789 if (! indexed_assign_conforms (iv.capacity (), 1, rhs_nr, rhs_nc)) |
|
3790 { |
|
3791 ::error ("A(matrix,int) = X: X must be a column vector with the"); |
|
3792 ::error ("same number of elements as matrix"); |
|
3793 return; |
|
3794 } |
|
3795 maybe_resize (iv.max (), j); |
|
3796 if (error_state) |
|
3797 return; |
|
3798 |
|
3799 do_matrix_assignment (rhs, iv, j); |
|
3800 } |
|
3801 break; |
|
3802 case complex_matrix_constant: |
|
3803 case matrix_constant: |
|
3804 { |
|
3805 Matrix mj = tmp_j.matrix_value (); |
|
3806 idx_vector jv (mj, user_pref.do_fortran_indexing, "column", |
|
3807 columns ()); |
|
3808 if (! jv) |
|
3809 return; |
|
3810 |
|
3811 if (! indexed_assign_conforms (iv.capacity (), jv.capacity (), |
|
3812 rhs_nr, rhs_nc)) |
|
3813 { |
|
3814 ::error ("A(r_mat,c_mat) = X: the number of rows in X must match"); |
|
3815 ::error ("the number of elements in r_mat and the number of"); |
|
3816 ::error ("columns in X must match the number of elements in c_mat"); |
|
3817 return; |
|
3818 } |
|
3819 maybe_resize (iv.max (), jv.max ()); |
|
3820 if (error_state) |
|
3821 return; |
|
3822 |
|
3823 do_matrix_assignment (rhs, iv, jv); |
|
3824 } |
|
3825 break; |
|
3826 case string_constant: |
|
3827 gripe_string_invalid (); |
|
3828 break; |
|
3829 case range_constant: |
|
3830 { |
|
3831 Range rj = tmp_j.range_value (); |
|
3832 if (! indexed_assign_conforms (iv.capacity (), rj.nelem (), |
|
3833 rhs_nr, rhs_nc)) |
|
3834 { |
|
3835 ::error ("A(matrix,range) = X: the number of rows in X must match"); |
|
3836 ::error ("the number of elements in matrix and the number of"); |
|
3837 ::error ("columns in X must match the number of elements in range"); |
|
3838 return; |
|
3839 } |
|
3840 |
|
3841 int nc = columns (); |
|
3842 if (nc == 2 && is_zero_one (rj) && rhs_nc == 1) |
|
3843 { |
|
3844 do_matrix_assignment (rhs, iv, 1); |
|
3845 } |
|
3846 else if (nc == 2 && is_one_zero (rj) && rhs_nc == 1) |
|
3847 { |
|
3848 do_matrix_assignment (rhs, iv, 0); |
|
3849 } |
|
3850 else |
|
3851 { |
|
3852 if (index_check (rj, "column") < 0) |
|
3853 return; |
|
3854 maybe_resize (iv.max (), tree_to_mat_idx (rj.max ())); |
|
3855 if (error_state) |
|
3856 return; |
|
3857 |
|
3858 do_matrix_assignment (rhs, iv, rj); |
|
3859 } |
|
3860 } |
|
3861 break; |
|
3862 case magic_colon: |
|
3863 { |
|
3864 int nc = columns (); |
|
3865 int new_nc = nc; |
|
3866 if (nc == 0) |
|
3867 new_nc = rhs_nc; |
|
3868 |
|
3869 if (indexed_assign_conforms (iv.capacity (), new_nc, |
|
3870 rhs_nr, rhs_nc)) |
|
3871 { |
|
3872 maybe_resize (iv.max (), new_nc-1); |
|
3873 if (error_state) |
|
3874 return; |
|
3875 } |
|
3876 else if (rhs_nr == 0 && rhs_nc == 0) |
|
3877 { |
|
3878 if (iv.max () >= rows ()) |
|
3879 { |
|
3880 ::error ("A(matrix,:) = []: row index out of range"); |
|
3881 return; |
|
3882 } |
|
3883 } |
|
3884 else |
|
3885 { |
|
3886 ::error ("A(matrix,:) = X: the number of rows in X must match the"); |
|
3887 ::error ("number of elements in matrix, and the number of columns"); |
|
3888 ::error ("in X must match the number of columns in A"); |
|
3889 return; |
|
3890 } |
|
3891 |
|
3892 do_matrix_assignment (rhs, iv, magic_colon); |
|
3893 } |
|
3894 break; |
|
3895 default: |
|
3896 panic_impossible (); |
|
3897 break; |
|
3898 } |
|
3899 } |
|
3900 |
|
3901 /* MA3 */ |
|
3902 void |
|
3903 tree_constant_rep::do_matrix_assignment (tree_constant& rhs, |
|
3904 Range& ri, tree_constant& j_arg) |
|
3905 { |
|
3906 tree_constant tmp_j = j_arg.make_numeric_or_range_or_magic (); |
|
3907 |
|
3908 tree_constant_rep::constant_type jtype = tmp_j.const_type (); |
|
3909 |
|
3910 int rhs_nr = rhs.rows (); |
|
3911 int rhs_nc = rhs.columns (); |
|
3912 |
|
3913 switch (jtype) |
|
3914 { |
|
3915 case complex_scalar_constant: |
|
3916 case scalar_constant: |
|
3917 { |
|
3918 int j = tree_to_mat_idx (tmp_j.double_value ()); |
|
3919 if (index_check (j, "column") < 0) |
|
3920 return; |
|
3921 if (! indexed_assign_conforms (ri.nelem (), 1, rhs_nr, rhs_nc)) |
|
3922 { |
|
3923 ::error ("A(range,int) = X: X must be a column vector with the"); |
|
3924 ::error ("same number of elements as range"); |
|
3925 return; |
|
3926 } |
|
3927 maybe_resize (tree_to_mat_idx (ri.max ()), j); |
|
3928 if (error_state) |
|
3929 return; |
|
3930 |
|
3931 do_matrix_assignment (rhs, ri, j); |
|
3932 } |
|
3933 break; |
|
3934 case complex_matrix_constant: |
|
3935 case matrix_constant: |
|
3936 { |
|
3937 Matrix mj = tmp_j.matrix_value (); |
|
3938 idx_vector jv (mj, user_pref.do_fortran_indexing, "column", |
|
3939 columns ()); |
|
3940 if (! jv) |
|
3941 return; |
|
3942 |
|
3943 if (! indexed_assign_conforms (ri.nelem (), jv.capacity (), |
|
3944 rhs_nr, rhs_nc)) |
|
3945 { |
|
3946 ::error ("A(range,matrix) = X: the number of rows in X must match"); |
|
3947 ::error ("the number of elements in range and the number of"); |
|
3948 ::error ("columns in X must match the number of elements in matrix"); |
|
3949 return; |
|
3950 } |
|
3951 maybe_resize (tree_to_mat_idx (ri.max ()), jv.max ()); |
|
3952 if (error_state) |
|
3953 return; |
|
3954 |
|
3955 do_matrix_assignment (rhs, ri, jv); |
|
3956 } |
|
3957 break; |
|
3958 case string_constant: |
|
3959 gripe_string_invalid (); |
|
3960 break; |
|
3961 case range_constant: |
|
3962 { |
|
3963 Range rj = tmp_j.range_value (); |
|
3964 if (! indexed_assign_conforms (ri.nelem (), rj.nelem (), |
|
3965 rhs_nr, rhs_nc)) |
|
3966 { |
|
3967 ::error ("A(r_range,c_range) = X: the number of rows in X must"); |
|
3968 ::error ("match the number of elements in r_range and the number"); |
|
3969 ::error ("of columns in X must match the number of elements in"); |
|
3970 ::error ("c_range"); |
|
3971 return; |
|
3972 } |
|
3973 |
|
3974 int nc = columns (); |
|
3975 if (nc == 2 && is_zero_one (rj) && rhs_nc == 1) |
|
3976 { |
|
3977 do_matrix_assignment (rhs, ri, 1); |
|
3978 } |
|
3979 else if (nc == 2 && is_one_zero (rj) && rhs_nc == 1) |
|
3980 { |
|
3981 do_matrix_assignment (rhs, ri, 0); |
|
3982 } |
|
3983 else |
|
3984 { |
|
3985 if (index_check (rj, "column") < 0) |
|
3986 return; |
|
3987 |
|
3988 maybe_resize (tree_to_mat_idx (ri.max ()), |
|
3989 tree_to_mat_idx (rj.max ())); |
|
3990 |
|
3991 if (error_state) |
|
3992 return; |
|
3993 |
|
3994 do_matrix_assignment (rhs, ri, rj); |
|
3995 } |
|
3996 } |
|
3997 break; |
|
3998 case magic_colon: |
|
3999 { |
|
4000 int nc = columns (); |
|
4001 int new_nc = nc; |
|
4002 if (nc == 0) |
|
4003 new_nc = rhs_nc; |
|
4004 |
|
4005 if (indexed_assign_conforms (ri.nelem (), new_nc, rhs_nr, rhs_nc)) |
|
4006 { |
|
4007 maybe_resize (tree_to_mat_idx (ri.max ()), new_nc-1); |
|
4008 if (error_state) |
|
4009 return; |
|
4010 } |
|
4011 else if (rhs_nr == 0 && rhs_nc == 0) |
|
4012 { |
|
4013 int b = tree_to_mat_idx (ri.min ()); |
|
4014 int l = tree_to_mat_idx (ri.max ()); |
|
4015 if (b < 0 || l >= rows ()) |
|
4016 { |
|
4017 ::error ("A(range,:) = []: row index out of range"); |
|
4018 return; |
|
4019 } |
|
4020 } |
|
4021 else |
|
4022 { |
|
4023 ::error ("A(range,:) = X: the number of rows in X must match the"); |
|
4024 ::error ("number of elements in range, and the number of columns"); |
|
4025 ::error ("in X must match the number of columns in A"); |
|
4026 return; |
|
4027 } |
|
4028 |
|
4029 do_matrix_assignment (rhs, ri, magic_colon); |
|
4030 } |
|
4031 break; |
|
4032 default: |
|
4033 panic_impossible (); |
|
4034 break; |
|
4035 } |
|
4036 } |
|
4037 |
|
4038 /* MA4 */ |
|
4039 void |
|
4040 tree_constant_rep::do_matrix_assignment (tree_constant& rhs, |
|
4041 tree_constant_rep::constant_type i, |
|
4042 tree_constant& j_arg) |
|
4043 { |
|
4044 tree_constant tmp_j = j_arg.make_numeric_or_range_or_magic (); |
|
4045 |
|
4046 tree_constant_rep::constant_type jtype = tmp_j.const_type (); |
|
4047 |
|
4048 int rhs_nr = rhs.rows (); |
|
4049 int rhs_nc = rhs.columns (); |
|
4050 |
|
4051 switch (jtype) |
|
4052 { |
|
4053 case complex_scalar_constant: |
|
4054 case scalar_constant: |
|
4055 { |
|
4056 int j = tree_to_mat_idx (tmp_j.double_value ()); |
|
4057 if (index_check (j, "column") < 0) |
|
4058 return; |
|
4059 int nr = rows (); |
|
4060 int nc = columns (); |
|
4061 if (nr == 0 && nc == 0 && rhs_nc == 1) |
|
4062 { |
|
4063 if (rhs.is_complex_type ()) |
|
4064 { |
|
4065 complex_matrix = new ComplexMatrix (); |
|
4066 type_tag = complex_matrix_constant; |
|
4067 } |
|
4068 else |
|
4069 { |
|
4070 matrix = new Matrix (); |
|
4071 type_tag = matrix_constant; |
|
4072 } |
|
4073 maybe_resize (rhs_nr-1, j); |
|
4074 if (error_state) |
|
4075 return; |
|
4076 } |
|
4077 else if (indexed_assign_conforms (nr, 1, rhs_nr, rhs_nc)) |
|
4078 { |
|
4079 maybe_resize (nr-1, j); |
|
4080 if (error_state) |
|
4081 return; |
|
4082 } |
|
4083 else if (rhs_nr == 0 && rhs_nc == 0) |
|
4084 { |
|
4085 if (j < 0 || j >= nc) |
|
4086 { |
|
4087 ::error ("A(:,int) = []: column index out of range"); |
|
4088 return; |
|
4089 } |
|
4090 } |
|
4091 else |
|
4092 { |
|
4093 ::error ("A(:,int) = X: X must be a column vector with the same"); |
|
4094 ::error ("number of rows as A"); |
|
4095 return; |
|
4096 } |
|
4097 |
|
4098 do_matrix_assignment (rhs, magic_colon, j); |
|
4099 } |
|
4100 break; |
|
4101 case complex_matrix_constant: |
|
4102 case matrix_constant: |
|
4103 { |
|
4104 Matrix mj = tmp_j.matrix_value (); |
|
4105 idx_vector jv (mj, user_pref.do_fortran_indexing, "column", |
|
4106 columns ()); |
|
4107 if (! jv) |
|
4108 return; |
|
4109 |
|
4110 int nr = rows (); |
|
4111 int new_nr = nr; |
|
4112 if (nr == 0) |
|
4113 new_nr = rhs_nr; |
|
4114 |
|
4115 if (indexed_assign_conforms (new_nr, jv.capacity (), |
|
4116 rhs_nr, rhs_nc)) |
|
4117 { |
|
4118 maybe_resize (new_nr-1, jv.max ()); |
|
4119 if (error_state) |
|
4120 return; |
|
4121 } |
|
4122 else if (rhs_nr == 0 && rhs_nc == 0) |
|
4123 { |
|
4124 if (jv.max () >= columns ()) |
|
4125 { |
|
4126 ::error ("A(:,matrix) = []: column index out of range"); |
|
4127 return; |
|
4128 } |
|
4129 } |
|
4130 else |
|
4131 { |
|
4132 ::error ("A(:,matrix) = X: the number of rows in X must match the"); |
|
4133 ::error ("number of rows in A, and the number of columns in X must"); |
|
4134 ::error ("match the number of elements in matrix"); |
|
4135 return; |
|
4136 } |
|
4137 |
|
4138 do_matrix_assignment (rhs, magic_colon, jv); |
|
4139 } |
|
4140 break; |
|
4141 case string_constant: |
|
4142 gripe_string_invalid (); |
|
4143 break; |
|
4144 case range_constant: |
|
4145 { |
|
4146 Range rj = tmp_j.range_value (); |
|
4147 int nr = rows (); |
|
4148 int new_nr = nr; |
|
4149 if (nr == 0) |
|
4150 new_nr = rhs_nr; |
|
4151 |
|
4152 if (indexed_assign_conforms (new_nr, rj.nelem (), rhs_nr, rhs_nc)) |
|
4153 { |
|
4154 int nc = columns (); |
|
4155 if (nc == 2 && is_zero_one (rj) && rhs_nc == 1) |
|
4156 { |
|
4157 do_matrix_assignment (rhs, magic_colon, 1); |
|
4158 } |
|
4159 else if (nc == 2 && is_one_zero (rj) && rhs_nc == 1) |
|
4160 { |
|
4161 do_matrix_assignment (rhs, magic_colon, 0); |
|
4162 } |
|
4163 else |
|
4164 { |
|
4165 if (index_check (rj, "column") < 0) |
|
4166 return; |
|
4167 maybe_resize (new_nr-1, tree_to_mat_idx (rj.max ())); |
|
4168 if (error_state) |
|
4169 return; |
|
4170 } |
|
4171 } |
|
4172 else if (rhs_nr == 0 && rhs_nc == 0) |
|
4173 { |
|
4174 int b = tree_to_mat_idx (rj.min ()); |
|
4175 int l = tree_to_mat_idx (rj.max ()); |
|
4176 if (b < 0 || l >= columns ()) |
|
4177 { |
|
4178 ::error ("A(:,range) = []: column index out of range"); |
|
4179 return; |
|
4180 } |
|
4181 } |
|
4182 else |
|
4183 { |
|
4184 ::error ("A(:,range) = X: the number of rows in X must match the"); |
|
4185 ::error ("number of rows in A, and the number of columns in X"); |
|
4186 ::error ("must match the number of elements in range"); |
|
4187 return; |
|
4188 } |
|
4189 |
|
4190 do_matrix_assignment (rhs, magic_colon, rj); |
|
4191 } |
|
4192 break; |
|
4193 case magic_colon: |
|
4194 // a(:,:) = foo is equivalent to a = foo. |
|
4195 do_matrix_assignment (rhs, magic_colon, magic_colon); |
|
4196 break; |
|
4197 default: |
|
4198 panic_impossible (); |
|
4199 break; |
|
4200 } |
|
4201 } |
|
4202 |
|
4203 /* |
|
4204 * Functions that actually handle assignment to a matrix using two |
|
4205 * index values. |
|
4206 * |
|
4207 * idx2 |
|
4208 * +---+---+----+----+ |
|
4209 * idx1 | i | v | r | c | |
|
4210 * ---------+---+---+----+----+ |
|
4211 * integer | 1 | 5 | 9 | 13 | |
|
4212 * ---------+---+---+----+----+ |
|
4213 * vector | 2 | 6 | 10 | 14 | |
|
4214 * ---------+---+---+----+----+ |
|
4215 * range | 3 | 7 | 11 | 15 | |
|
4216 * ---------+---+---+----+----+ |
|
4217 * colon | 4 | 8 | 12 | 16 | |
|
4218 * ---------+---+---+----+----+ |
|
4219 */ |
|
4220 |
|
4221 /* 1 */ |
|
4222 void |
|
4223 tree_constant_rep::do_matrix_assignment (tree_constant& rhs, int i, int j) |
|
4224 { |
|
4225 REP_ELEM_ASSIGN (i, j, rhs.double_value (), rhs.complex_value (), |
|
4226 rhs.is_real_type ()); |
|
4227 } |
|
4228 |
|
4229 /* 2 */ |
|
4230 void |
|
4231 tree_constant_rep::do_matrix_assignment (tree_constant& rhs, int i, |
|
4232 idx_vector& jv) |
|
4233 { |
|
4234 REP_RHS_MATRIX (rhs, rhs_m, rhs_cm, rhs_nr, rhs_nc); |
|
4235 |
|
4236 for (int j = 0; j < jv.capacity (); j++) |
|
4237 REP_ELEM_ASSIGN (i, jv.elem (j), rhs_m.elem (0, j), |
|
4238 rhs_cm.elem (0, j), rhs.is_real_type ()); |
|
4239 } |
|
4240 |
|
4241 /* 3 */ |
|
4242 void |
|
4243 tree_constant_rep::do_matrix_assignment (tree_constant& rhs, int i, Range& rj) |
|
4244 { |
|
4245 REP_RHS_MATRIX (rhs, rhs_m, rhs_cm, rhs_nr, rhs_nc); |
|
4246 |
|
4247 double b = rj.base (); |
|
4248 double increment = rj.inc (); |
|
4249 |
|
4250 for (int j = 0; j < rj.nelem (); j++) |
|
4251 { |
|
4252 double tmp = b + j * increment; |
|
4253 int col = tree_to_mat_idx (tmp); |
|
4254 REP_ELEM_ASSIGN (i, col, rhs_m.elem (0, j), rhs_cm.elem (0, j), |
|
4255 rhs.is_real_type ()); |
|
4256 } |
|
4257 } |
|
4258 |
|
4259 /* 4 */ |
|
4260 void |
|
4261 tree_constant_rep::do_matrix_assignment (tree_constant& rhs, int i, |
|
4262 tree_constant_rep::constant_type mcj) |
|
4263 { |
|
4264 assert (mcj == magic_colon); |
|
4265 |
|
4266 int nc = columns (); |
|
4267 |
|
4268 if (rhs.is_zero_by_zero ()) |
|
4269 { |
|
4270 delete_row (i); |
|
4271 } |
|
4272 else if (rhs.is_matrix_type ()) |
|
4273 { |
|
4274 REP_RHS_MATRIX (rhs, rhs_m, rhs_cm, rhs_nr, rhs_nc); |
|
4275 |
|
4276 for (int j = 0; j < nc; j++) |
|
4277 REP_ELEM_ASSIGN (i, j, rhs_m.elem (0, j), rhs_cm.elem (0, j), |
|
4278 rhs.is_real_type ()); |
|
4279 } |
|
4280 else if (rhs.is_scalar_type () && nc == 1) |
|
4281 { |
|
4282 REP_ELEM_ASSIGN (i, 0, rhs.double_value (), |
|
4283 rhs.complex_value (), rhs.is_real_type ()); |
|
4284 } |
|
4285 else |
|
4286 panic_impossible (); |
|
4287 } |
|
4288 |
|
4289 /* 5 */ |
|
4290 void |
|
4291 tree_constant_rep::do_matrix_assignment (tree_constant& rhs, |
|
4292 idx_vector& iv, int j) |
|
4293 { |
|
4294 REP_RHS_MATRIX (rhs, rhs_m, rhs_cm, rhs_nr, rhs_nc); |
|
4295 |
|
4296 for (int i = 0; i < iv.capacity (); i++) |
|
4297 { |
|
4298 int row = iv.elem (i); |
|
4299 REP_ELEM_ASSIGN (row, j, rhs_m.elem (i, 0), |
|
4300 rhs_cm.elem (i, 0), rhs.is_real_type ()); |
|
4301 } |
|
4302 } |
|
4303 |
|
4304 /* 6 */ |
|
4305 void |
|
4306 tree_constant_rep::do_matrix_assignment (tree_constant& rhs, |
|
4307 idx_vector& iv, idx_vector& jv) |
|
4308 { |
|
4309 REP_RHS_MATRIX (rhs, rhs_m, rhs_cm, rhs_nr, rhs_nc); |
|
4310 |
|
4311 for (int i = 0; i < iv.capacity (); i++) |
|
4312 { |
|
4313 int row = iv.elem (i); |
|
4314 for (int j = 0; j < jv.capacity (); j++) |
|
4315 { |
|
4316 int col = jv.elem (j); |
|
4317 REP_ELEM_ASSIGN (row, col, rhs_m.elem (i, j), |
|
4318 rhs_cm.elem (i, j), rhs.is_real_type ()); |
|
4319 } |
|
4320 } |
|
4321 } |
|
4322 |
|
4323 /* 7 */ |
|
4324 void |
|
4325 tree_constant_rep::do_matrix_assignment (tree_constant& rhs, |
|
4326 idx_vector& iv, Range& rj) |
|
4327 { |
|
4328 REP_RHS_MATRIX (rhs, rhs_m, rhs_cm, rhs_nr, rhs_nc); |
|
4329 |
|
4330 double b = rj.base (); |
|
4331 double increment = rj.inc (); |
|
4332 |
|
4333 for (int i = 0; i < iv.capacity (); i++) |
|
4334 { |
|
4335 int row = iv.elem (i); |
|
4336 for (int j = 0; j < rj.nelem (); j++) |
|
4337 { |
|
4338 double tmp = b + j * increment; |
|
4339 int col = tree_to_mat_idx (tmp); |
|
4340 REP_ELEM_ASSIGN (row, col, rhs_m.elem (i, j), |
|
4341 rhs_cm.elem (i, j), rhs.is_real_type ()); |
|
4342 } |
|
4343 } |
|
4344 } |
|
4345 |
|
4346 /* 8 */ |
|
4347 void |
|
4348 tree_constant_rep::do_matrix_assignment (tree_constant& rhs, idx_vector& iv, |
|
4349 tree_constant_rep::constant_type mcj) |
|
4350 { |
|
4351 assert (mcj == magic_colon); |
|
4352 |
|
4353 if (rhs.is_zero_by_zero ()) |
|
4354 { |
|
4355 delete_rows (iv); |
|
4356 } |
|
4357 else |
|
4358 { |
|
4359 REP_RHS_MATRIX (rhs, rhs_m, rhs_cm, rhs_nr, rhs_nc); |
|
4360 |
|
4361 int nc = columns (); |
|
4362 |
|
4363 for (int j = 0; j < nc; j++) |
|
4364 { |
|
4365 for (int i = 0; i < iv.capacity (); i++) |
|
4366 { |
|
4367 int row = iv.elem (i); |
|
4368 REP_ELEM_ASSIGN (row, j, rhs_m.elem (i, j), |
|
4369 rhs_cm.elem (i, j), rhs.is_real_type ()); |
|
4370 } |
|
4371 } |
|
4372 } |
|
4373 } |
|
4374 |
|
4375 /* 9 */ |
|
4376 void |
|
4377 tree_constant_rep::do_matrix_assignment (tree_constant& rhs, Range& ri, int j) |
|
4378 { |
|
4379 REP_RHS_MATRIX (rhs, rhs_m, rhs_cm, rhs_nr, rhs_nc); |
|
4380 |
|
4381 double b = ri.base (); |
|
4382 double increment = ri.inc (); |
|
4383 |
|
4384 for (int i = 0; i < ri.nelem (); i++) |
|
4385 { |
|
4386 double tmp = b + i * increment; |
|
4387 int row = tree_to_mat_idx (tmp); |
|
4388 REP_ELEM_ASSIGN (row, j, rhs_m.elem (i, 0), |
|
4389 rhs_cm.elem (i, 0), rhs.is_real_type ()); |
|
4390 } |
|
4391 } |
|
4392 |
|
4393 /* 10 */ |
|
4394 void |
|
4395 tree_constant_rep::do_matrix_assignment (tree_constant& rhs, Range& ri, |
|
4396 idx_vector& jv) |
|
4397 { |
|
4398 REP_RHS_MATRIX (rhs, rhs_m, rhs_cm, rhs_nr, rhs_nc); |
|
4399 |
|
4400 double b = ri.base (); |
|
4401 double increment = ri.inc (); |
|
4402 |
|
4403 for (int j = 0; j < jv.capacity (); j++) |
|
4404 { |
|
4405 int col = jv.elem (j); |
|
4406 for (int i = 0; i < ri.nelem (); i++) |
|
4407 { |
|
4408 double tmp = b + i * increment; |
|
4409 int row = tree_to_mat_idx (tmp); |
|
4410 REP_ELEM_ASSIGN (row, col, rhs_m.elem (i, j), |
|
4411 rhs_m.elem (i, j), rhs.is_real_type ()); |
|
4412 } |
|
4413 } |
|
4414 } |
|
4415 |
|
4416 /* 11 */ |
|
4417 void |
|
4418 tree_constant_rep::do_matrix_assignment (tree_constant& rhs, Range& ri, |
|
4419 Range& rj) |
|
4420 { |
|
4421 double ib = ri.base (); |
|
4422 double iinc = ri.inc (); |
|
4423 double jb = rj.base (); |
|
4424 double jinc = rj.inc (); |
|
4425 |
|
4426 REP_RHS_MATRIX (rhs, rhs_m, rhs_cm, rhs_nr, rhs_nc); |
|
4427 |
|
4428 for (int i = 0; i < ri.nelem (); i++) |
|
4429 { |
|
4430 double itmp = ib + i * iinc; |
|
4431 int row = tree_to_mat_idx (itmp); |
|
4432 for (int j = 0; j < rj.nelem (); j++) |
|
4433 { |
|
4434 double jtmp = jb + j * jinc; |
|
4435 int col = tree_to_mat_idx (jtmp); |
|
4436 REP_ELEM_ASSIGN (row, col, rhs_m.elem (i, j), |
|
4437 rhs_cm.elem (i, j), rhs.is_real_type ()); |
|
4438 } |
|
4439 } |
|
4440 } |
|
4441 |
|
4442 /* 12 */ |
|
4443 void |
|
4444 tree_constant_rep::do_matrix_assignment (tree_constant& rhs, Range& ri, |
|
4445 tree_constant_rep::constant_type mcj) |
|
4446 { |
|
4447 assert (mcj == magic_colon); |
|
4448 |
|
4449 if (rhs.is_zero_by_zero ()) |
|
4450 { |
|
4451 delete_rows (ri); |
|
4452 } |
|
4453 else |
|
4454 { |
|
4455 REP_RHS_MATRIX (rhs, rhs_m, rhs_cm, rhs_nr, rhs_nc); |
|
4456 |
|
4457 double ib = ri.base (); |
|
4458 double iinc = ri.inc (); |
|
4459 |
|
4460 int nc = columns (); |
|
4461 |
|
4462 for (int i = 0; i < ri.nelem (); i++) |
|
4463 { |
|
4464 double itmp = ib + i * iinc; |
|
4465 int row = tree_to_mat_idx (itmp); |
|
4466 for (int j = 0; j < nc; j++) |
|
4467 REP_ELEM_ASSIGN (row, j, rhs_m.elem (i, j), |
|
4468 rhs_cm.elem (i, j), rhs.is_real_type ()); |
|
4469 } |
|
4470 } |
|
4471 } |
|
4472 |
|
4473 /* 13 */ |
|
4474 void |
|
4475 tree_constant_rep::do_matrix_assignment (tree_constant& rhs, |
|
4476 tree_constant_rep::constant_type mci, |
|
4477 int j) |
|
4478 { |
|
4479 assert (mci == magic_colon); |
|
4480 |
|
4481 int nr = rows (); |
|
4482 |
|
4483 if (rhs.is_zero_by_zero ()) |
|
4484 { |
|
4485 delete_column (j); |
|
4486 } |
|
4487 else if (rhs.is_matrix_type ()) |
|
4488 { |
|
4489 REP_RHS_MATRIX (rhs, rhs_m, rhs_cm, rhs_nr, rhs_nc); |
|
4490 |
|
4491 for (int i = 0; i < nr; i++) |
|
4492 REP_ELEM_ASSIGN (i, j, rhs_m.elem (i, 0), |
|
4493 rhs_cm.elem (i, 0), rhs.is_real_type ()); |
|
4494 } |
|
4495 else if (rhs.is_scalar_type () && nr == 1) |
|
4496 { |
|
4497 REP_ELEM_ASSIGN (0, j, rhs.double_value (), |
|
4498 rhs.complex_value (), rhs.is_real_type ()); |
|
4499 } |
|
4500 else |
|
4501 panic_impossible (); |
|
4502 } |
|
4503 |
|
4504 /* 14 */ |
|
4505 void |
|
4506 tree_constant_rep::do_matrix_assignment (tree_constant& rhs, |
|
4507 tree_constant_rep::constant_type mci, |
|
4508 idx_vector& jv) |
|
4509 { |
|
4510 assert (mci == magic_colon); |
|
4511 |
|
4512 if (rhs.is_zero_by_zero ()) |
|
4513 { |
|
4514 delete_columns (jv); |
|
4515 } |
|
4516 else |
|
4517 { |
|
4518 REP_RHS_MATRIX (rhs, rhs_m, rhs_cm, rhs_nr, rhs_nc); |
|
4519 |
|
4520 int nr = rows (); |
|
4521 |
|
4522 for (int i = 0; i < nr; i++) |
|
4523 { |
|
4524 for (int j = 0; j < jv.capacity (); j++) |
|
4525 { |
|
4526 int col = jv.elem (j); |
|
4527 REP_ELEM_ASSIGN (i, col, rhs_m.elem (i, j), |
|
4528 rhs_cm.elem (i, j), rhs.is_real_type ()); |
|
4529 } |
|
4530 } |
|
4531 } |
|
4532 } |
|
4533 |
|
4534 /* 15 */ |
|
4535 void |
|
4536 tree_constant_rep::do_matrix_assignment (tree_constant& rhs, |
|
4537 tree_constant_rep::constant_type mci, |
|
4538 Range& rj) |
|
4539 { |
|
4540 assert (mci == magic_colon); |
|
4541 |
|
4542 if (rhs.is_zero_by_zero ()) |
|
4543 { |
|
4544 delete_columns (rj); |
|
4545 } |
|
4546 else |
|
4547 { |
|
4548 REP_RHS_MATRIX (rhs, rhs_m, rhs_cm, rhs_nr, rhs_nc); |
|
4549 |
|
4550 int nr = rows (); |
|
4551 |
|
4552 double jb = rj.base (); |
|
4553 double jinc = rj.inc (); |
|
4554 |
|
4555 for (int j = 0; j < rj.nelem (); j++) |
|
4556 { |
|
4557 double jtmp = jb + j * jinc; |
|
4558 int col = tree_to_mat_idx (jtmp); |
|
4559 for (int i = 0; i < nr; i++) |
|
4560 { |
|
4561 REP_ELEM_ASSIGN (i, col, rhs_m.elem (i, j), |
|
4562 rhs_cm.elem (i, j), rhs.is_real_type ()); |
|
4563 } |
|
4564 } |
|
4565 } |
|
4566 } |
|
4567 |
|
4568 /* 16 */ |
|
4569 void |
|
4570 tree_constant_rep::do_matrix_assignment (tree_constant& rhs, |
|
4571 tree_constant_rep::constant_type mci, |
|
4572 tree_constant_rep::constant_type mcj) |
|
4573 { |
|
4574 assert (mci == magic_colon && mcj == magic_colon); |
|
4575 |
|
4576 switch (type_tag) |
|
4577 { |
|
4578 case scalar_constant: |
|
4579 break; |
|
4580 case matrix_constant: |
|
4581 delete matrix; |
|
4582 break; |
|
4583 case complex_scalar_constant: |
|
4584 delete complex_scalar; |
|
4585 break; |
|
4586 case complex_matrix_constant: |
|
4587 delete complex_matrix; |
|
4588 break; |
|
4589 case string_constant: |
|
4590 delete [] string; |
|
4591 break; |
|
4592 case range_constant: |
|
4593 delete range; |
|
4594 break; |
|
4595 case magic_colon: |
|
4596 default: |
|
4597 panic_impossible (); |
|
4598 break; |
|
4599 } |
|
4600 |
|
4601 type_tag = rhs.const_type (); |
|
4602 |
|
4603 switch (type_tag) |
|
4604 { |
|
4605 case scalar_constant: |
|
4606 scalar = rhs.double_value (); |
|
4607 break; |
|
4608 case matrix_constant: |
|
4609 matrix = new Matrix (rhs.matrix_value ()); |
|
4610 break; |
|
4611 case string_constant: |
|
4612 string = strsave (rhs.string_value ()); |
|
4613 break; |
|
4614 case complex_matrix_constant: |
|
4615 complex_matrix = new ComplexMatrix (rhs.complex_matrix_value ()); |
|
4616 break; |
|
4617 case complex_scalar_constant: |
|
4618 complex_scalar = new Complex (rhs.complex_value ()); |
|
4619 break; |
|
4620 case range_constant: |
|
4621 range = new Range (rhs.range_value ()); |
|
4622 break; |
|
4623 case magic_colon: |
|
4624 default: |
|
4625 panic_impossible (); |
|
4626 break; |
|
4627 } |
|
4628 } |
|
4629 |
|
4630 /* |
|
4631 * Functions for deleting rows or columns of a matrix. These are used |
|
4632 * to handle statements like |
|
4633 * |
|
4634 * M (i, j) = [] |
|
4635 */ |
|
4636 void |
|
4637 tree_constant_rep::delete_row (int idx) |
|
4638 { |
|
4639 if (type_tag == matrix_constant) |
|
4640 { |
|
4641 int nr = matrix->rows (); |
|
4642 int nc = matrix->columns (); |
|
4643 Matrix *new_matrix = new Matrix (nr-1, nc); |
|
4644 int ii = 0; |
|
4645 for (int i = 0; i < nr; i++) |
|
4646 { |
|
4647 if (i != idx) |
|
4648 { |
|
4649 for (int j = 0; j < nc; j++) |
|
4650 new_matrix->elem (ii, j) = matrix->elem (i, j); |
|
4651 ii++; |
|
4652 } |
|
4653 } |
|
4654 delete matrix; |
|
4655 matrix = new_matrix; |
|
4656 } |
|
4657 else if (type_tag == complex_matrix_constant) |
|
4658 { |
|
4659 int nr = complex_matrix->rows (); |
|
4660 int nc = complex_matrix->columns (); |
|
4661 ComplexMatrix *new_matrix = new ComplexMatrix (nr-1, nc); |
|
4662 int ii = 0; |
|
4663 for (int i = 0; i < nr; i++) |
|
4664 { |
|
4665 if (i != idx) |
|
4666 { |
|
4667 for (int j = 0; j < nc; j++) |
|
4668 new_matrix->elem (ii, j) = complex_matrix->elem (i, j); |
|
4669 ii++; |
|
4670 } |
|
4671 } |
|
4672 delete complex_matrix; |
|
4673 complex_matrix = new_matrix; |
|
4674 } |
|
4675 else |
|
4676 panic_impossible (); |
|
4677 } |
|
4678 |
|
4679 void |
|
4680 tree_constant_rep::delete_rows (idx_vector& iv) |
|
4681 { |
|
4682 iv.sort_uniq (); |
|
4683 int num_to_delete = iv.length (); |
|
4684 |
|
4685 int nr = rows (); |
|
4686 int nc = columns (); |
|
4687 |
|
4688 // If deleting all rows of a column vector, make result 0x0. |
|
4689 if (nc == 1 && num_to_delete == nr) |
|
4690 nc = 0; |
|
4691 |
|
4692 if (type_tag == matrix_constant) |
|
4693 { |
|
4694 Matrix *new_matrix = new Matrix (nr-num_to_delete, nc); |
|
4695 if (nr > num_to_delete) |
|
4696 { |
|
4697 int ii = 0; |
|
4698 int idx = 0; |
|
4699 for (int i = 0; i < nr; i++) |
|
4700 { |
|
4701 if (i == iv.elem (idx)) |
|
4702 idx++; |
|
4703 else |
|
4704 { |
|
4705 for (int j = 0; j < nc; j++) |
|
4706 new_matrix->elem (ii, j) = matrix->elem (i, j); |
|
4707 ii++; |
|
4708 } |
|
4709 } |
|
4710 } |
|
4711 delete matrix; |
|
4712 matrix = new_matrix; |
|
4713 } |
|
4714 else if (type_tag == complex_matrix_constant) |
|
4715 { |
|
4716 ComplexMatrix *new_matrix = new ComplexMatrix (nr-num_to_delete, nc); |
|
4717 if (nr > num_to_delete) |
|
4718 { |
|
4719 int ii = 0; |
|
4720 int idx = 0; |
|
4721 for (int i = 0; i < nr; i++) |
|
4722 { |
|
4723 if (i == iv.elem (idx)) |
|
4724 idx++; |
|
4725 else |
|
4726 { |
|
4727 for (int j = 0; j < nc; j++) |
|
4728 new_matrix->elem (ii, j) = complex_matrix->elem (i, j); |
|
4729 ii++; |
|
4730 } |
|
4731 } |
|
4732 } |
|
4733 delete complex_matrix; |
|
4734 complex_matrix = new_matrix; |
|
4735 } |
|
4736 else |
|
4737 panic_impossible (); |
|
4738 } |
|
4739 |
|
4740 void |
|
4741 tree_constant_rep::delete_rows (Range& ri) |
|
4742 { |
|
4743 ri.sort (); |
|
4744 int num_to_delete = ri.nelem (); |
|
4745 |
|
4746 int nr = rows (); |
|
4747 int nc = columns (); |
|
4748 |
|
4749 // If deleting all rows of a column vector, make result 0x0. |
|
4750 if (nc == 1 && num_to_delete == nr) |
|
4751 nc = 0; |
|
4752 |
|
4753 double ib = ri.base (); |
|
4754 double iinc = ri.inc (); |
|
4755 |
|
4756 int max_idx = tree_to_mat_idx (ri.max ()); |
|
4757 |
|
4758 if (type_tag == matrix_constant) |
|
4759 { |
|
4760 Matrix *new_matrix = new Matrix (nr-num_to_delete, nc); |
|
4761 if (nr > num_to_delete) |
|
4762 { |
|
4763 int ii = 0; |
|
4764 int idx = 0; |
|
4765 for (int i = 0; i < nr; i++) |
|
4766 { |
|
4767 double itmp = ib + idx * iinc; |
|
4768 int row = tree_to_mat_idx (itmp); |
|
4769 |
|
4770 if (i == row && row <= max_idx) |
|
4771 idx++; |
|
4772 else |
|
4773 { |
|
4774 for (int j = 0; j < nc; j++) |
|
4775 new_matrix->elem (ii, j) = matrix->elem (i, j); |
|
4776 ii++; |
|
4777 } |
|
4778 } |
|
4779 } |
|
4780 delete matrix; |
|
4781 matrix = new_matrix; |
|
4782 } |
|
4783 else if (type_tag == complex_matrix_constant) |
|
4784 { |
|
4785 ComplexMatrix *new_matrix = new ComplexMatrix (nr-num_to_delete, nc); |
|
4786 if (nr > num_to_delete) |
|
4787 { |
|
4788 int ii = 0; |
|
4789 int idx = 0; |
|
4790 for (int i = 0; i < nr; i++) |
|
4791 { |
|
4792 double itmp = ib + idx * iinc; |
|
4793 int row = tree_to_mat_idx (itmp); |
|
4794 |
|
4795 if (i == row && row <= max_idx) |
|
4796 idx++; |
|
4797 else |
|
4798 { |
|
4799 for (int j = 0; j < nc; j++) |
|
4800 new_matrix->elem (ii, j) = complex_matrix->elem (i, j); |
|
4801 ii++; |
|
4802 } |
|
4803 } |
|
4804 } |
|
4805 delete complex_matrix; |
|
4806 complex_matrix = new_matrix; |
|
4807 } |
|
4808 else |
|
4809 panic_impossible (); |
|
4810 } |
|
4811 |
|
4812 void |
|
4813 tree_constant_rep::delete_column (int idx) |
|
4814 { |
|
4815 if (type_tag == matrix_constant) |
|
4816 { |
|
4817 int nr = matrix->rows (); |
|
4818 int nc = matrix->columns (); |
|
4819 Matrix *new_matrix = new Matrix (nr, nc-1); |
|
4820 int jj = 0; |
|
4821 for (int j = 0; j < nc; j++) |
|
4822 { |
|
4823 if (j != idx) |
|
4824 { |
|
4825 for (int i = 0; i < nr; i++) |
|
4826 new_matrix->elem (i, jj) = matrix->elem (i, j); |
|
4827 jj++; |
|
4828 } |
|
4829 } |
|
4830 delete matrix; |
|
4831 matrix = new_matrix; |
|
4832 } |
|
4833 else if (type_tag == complex_matrix_constant) |
|
4834 { |
|
4835 int nr = complex_matrix->rows (); |
|
4836 int nc = complex_matrix->columns (); |
|
4837 ComplexMatrix *new_matrix = new ComplexMatrix (nr, nc-1); |
|
4838 int jj = 0; |
|
4839 for (int j = 0; j < nc; j++) |
|
4840 { |
|
4841 if (j != idx) |
|
4842 { |
|
4843 for (int i = 0; i < nr; i++) |
|
4844 new_matrix->elem (i, jj) = complex_matrix->elem (i, j); |
|
4845 jj++; |
|
4846 } |
|
4847 } |
|
4848 delete complex_matrix; |
|
4849 complex_matrix = new_matrix; |
|
4850 } |
|
4851 else |
|
4852 panic_impossible (); |
|
4853 } |
|
4854 |
|
4855 void |
|
4856 tree_constant_rep::delete_columns (idx_vector& jv) |
|
4857 { |
|
4858 jv.sort_uniq (); |
|
4859 int num_to_delete = jv.length (); |
|
4860 |
|
4861 int nr = rows (); |
|
4862 int nc = columns (); |
|
4863 |
|
4864 // If deleting all columns of a row vector, make result 0x0. |
|
4865 if (nr == 1 && num_to_delete == nc) |
|
4866 nr = 0; |
|
4867 |
|
4868 if (type_tag == matrix_constant) |
|
4869 { |
|
4870 Matrix *new_matrix = new Matrix (nr, nc-num_to_delete); |
|
4871 if (nc > num_to_delete) |
|
4872 { |
|
4873 int jj = 0; |
|
4874 int idx = 0; |
|
4875 for (int j = 0; j < nc; j++) |
|
4876 { |
|
4877 if (j == jv.elem (idx)) |
|
4878 idx++; |
|
4879 else |
|
4880 { |
|
4881 for (int i = 0; i < nr; i++) |
|
4882 new_matrix->elem (i, jj) = matrix->elem (i, j); |
|
4883 jj++; |
|
4884 } |
|
4885 } |
|
4886 } |
|
4887 delete matrix; |
|
4888 matrix = new_matrix; |
|
4889 } |
|
4890 else if (type_tag == complex_matrix_constant) |
|
4891 { |
|
4892 ComplexMatrix *new_matrix = new ComplexMatrix (nr, nc-num_to_delete); |
|
4893 if (nc > num_to_delete) |
|
4894 { |
|
4895 int jj = 0; |
|
4896 int idx = 0; |
|
4897 for (int j = 0; j < nc; j++) |
|
4898 { |
|
4899 if (j == jv.elem (idx)) |
|
4900 idx++; |
|
4901 else |
|
4902 { |
|
4903 for (int i = 0; i < nr; i++) |
|
4904 new_matrix->elem (i, jj) = complex_matrix->elem (i, j); |
|
4905 jj++; |
|
4906 } |
|
4907 } |
|
4908 } |
|
4909 delete complex_matrix; |
|
4910 complex_matrix = new_matrix; |
|
4911 } |
|
4912 else |
|
4913 panic_impossible (); |
|
4914 } |
|
4915 |
|
4916 void |
|
4917 tree_constant_rep::delete_columns (Range& rj) |
|
4918 { |
|
4919 rj.sort (); |
|
4920 int num_to_delete = rj.nelem (); |
|
4921 |
|
4922 int nr = rows (); |
|
4923 int nc = columns (); |
|
4924 |
|
4925 // If deleting all columns of a row vector, make result 0x0. |
|
4926 if (nr == 1 && num_to_delete == nc) |
|
4927 nr = 0; |
|
4928 |
|
4929 double jb = rj.base (); |
|
4930 double jinc = rj.inc (); |
|
4931 |
|
4932 int max_idx = tree_to_mat_idx (rj.max ()); |
|
4933 |
|
4934 if (type_tag == matrix_constant) |
|
4935 { |
|
4936 Matrix *new_matrix = new Matrix (nr, nc-num_to_delete); |
|
4937 if (nc > num_to_delete) |
|
4938 { |
|
4939 int jj = 0; |
|
4940 int idx = 0; |
|
4941 for (int j = 0; j < nc; j++) |
|
4942 { |
|
4943 double jtmp = jb + idx * jinc; |
|
4944 int col = tree_to_mat_idx (jtmp); |
|
4945 |
|
4946 if (j == col && col <= max_idx) |
|
4947 idx++; |
|
4948 else |
|
4949 { |
|
4950 for (int i = 0; i < nr; i++) |
|
4951 new_matrix->elem (i, jj) = matrix->elem (i, j); |
|
4952 jj++; |
|
4953 } |
|
4954 } |
|
4955 } |
|
4956 delete matrix; |
|
4957 matrix = new_matrix; |
|
4958 } |
|
4959 else if (type_tag == complex_matrix_constant) |
|
4960 { |
|
4961 ComplexMatrix *new_matrix = new ComplexMatrix (nr, nc-num_to_delete); |
|
4962 if (nc > num_to_delete) |
|
4963 { |
|
4964 int jj = 0; |
|
4965 int idx = 0; |
|
4966 for (int j = 0; j < nc; j++) |
|
4967 { |
|
4968 double jtmp = jb + idx * jinc; |
|
4969 int col = tree_to_mat_idx (jtmp); |
|
4970 |
|
4971 if (j == col && col <= max_idx) |
|
4972 idx++; |
|
4973 else |
|
4974 { |
|
4975 for (int i = 0; i < nr; i++) |
|
4976 new_matrix->elem (i, jj) = complex_matrix->elem (i, j); |
|
4977 jj++; |
|
4978 } |
|
4979 } |
|
4980 } |
|
4981 delete complex_matrix; |
|
4982 complex_matrix = new_matrix; |
|
4983 } |
|
4984 else |
|
4985 panic_impossible (); |
|
4986 } |
|
4987 |
|
4988 /* |
|
4989 * Indexing functions. |
|
4990 */ |
|
4991 int |
|
4992 tree_constant_rep::valid_as_scalar_index (void) const |
|
4993 { |
|
4994 int valid = type_tag == magic_colon |
|
4995 || (type_tag == scalar_constant && NINT (scalar) == 1) |
|
4996 || (type_tag == range_constant |
|
4997 && range->nelem () == 1 && NINT (range->base ()) == 1); |
|
4998 |
|
4999 return valid; |
|
5000 } |
|
5001 |
|
5002 tree_constant |
500
|
5003 tree_constant_rep::do_scalar_index (const Octave_object& args, |
492
|
5004 int nargs) const |
|
5005 { |
|
5006 if (valid_scalar_indices (args, nargs)) |
|
5007 { |
|
5008 if (type_tag == scalar_constant) |
|
5009 return tree_constant (scalar); |
|
5010 else if (type_tag == complex_scalar_constant) |
|
5011 return tree_constant (*complex_scalar); |
|
5012 else |
|
5013 panic_impossible (); |
|
5014 } |
|
5015 else |
|
5016 { |
|
5017 int rows = 0; |
|
5018 int cols = 0; |
|
5019 |
|
5020 switch (nargs) |
|
5021 { |
|
5022 case 3: |
|
5023 { |
500
|
5024 if (args(2).is_matrix_type ()) |
492
|
5025 { |
500
|
5026 Matrix mj = args(2).matrix_value (); |
492
|
5027 |
|
5028 idx_vector j (mj, user_pref.do_fortran_indexing, ""); |
|
5029 if (! j) |
|
5030 return tree_constant (); |
|
5031 |
|
5032 int len = j.length (); |
|
5033 if (len == j.ones_count ()) |
|
5034 cols = len; |
|
5035 } |
500
|
5036 else if (args(2).const_type () == magic_colon |
|
5037 || (args(2).is_scalar_type () |
|
5038 && NINT (args(2).double_value ()) == 1)) |
492
|
5039 { |
|
5040 cols = 1; |
|
5041 } |
|
5042 else |
|
5043 break; |
|
5044 } |
|
5045 // Fall through... |
|
5046 case 2: |
|
5047 { |
500
|
5048 if (args(1).is_matrix_type ()) |
492
|
5049 { |
500
|
5050 Matrix mi = args(1).matrix_value (); |
492
|
5051 |
|
5052 idx_vector i (mi, user_pref.do_fortran_indexing, ""); |
|
5053 if (! i) |
|
5054 return tree_constant (); |
|
5055 |
|
5056 int len = i.length (); |
|
5057 if (len == i.ones_count ()) |
|
5058 rows = len; |
|
5059 } |
500
|
5060 else if (args(1).const_type () == magic_colon |
|
5061 || (args(1).is_scalar_type () |
|
5062 && NINT (args(1).double_value ()) == 1)) |
492
|
5063 { |
|
5064 rows = 1; |
|
5065 } |
500
|
5066 else if (args(1).is_scalar_type () |
|
5067 && NINT (args(1).double_value ()) == 0) |
492
|
5068 { |
|
5069 Matrix m (0, 0); |
|
5070 return tree_constant (m); |
|
5071 } |
|
5072 else |
|
5073 break; |
|
5074 |
|
5075 if (cols == 0) |
|
5076 { |
|
5077 if (user_pref.prefer_column_vectors) |
|
5078 cols = 1; |
|
5079 else |
|
5080 { |
|
5081 cols = rows; |
|
5082 rows = 1; |
|
5083 } |
|
5084 } |
|
5085 |
|
5086 if (type_tag == scalar_constant) |
|
5087 { |
|
5088 Matrix m (rows, cols, scalar); |
|
5089 return tree_constant (m); |
|
5090 } |
|
5091 else if (type_tag == complex_scalar_constant) |
|
5092 { |
|
5093 ComplexMatrix cm (rows, cols, *complex_scalar); |
|
5094 return tree_constant (cm); |
|
5095 } |
|
5096 else |
|
5097 panic_impossible (); |
|
5098 } |
|
5099 break; |
|
5100 default: |
|
5101 ::error ("illegal number of arguments for scalar type"); |
|
5102 return tree_constant (); |
|
5103 break; |
|
5104 } |
|
5105 } |
|
5106 |
|
5107 ::error ("index invalid or out of range for scalar type"); |
|
5108 return tree_constant (); |
|
5109 } |
|
5110 |
|
5111 tree_constant |
500
|
5112 tree_constant_rep::do_matrix_index (const Octave_object& args, |
492
|
5113 int nargin) const |
|
5114 { |
|
5115 tree_constant retval; |
|
5116 |
|
5117 switch (nargin) |
|
5118 { |
|
5119 case 2: |
500
|
5120 if (args.length () <= 0) |
492
|
5121 ::error ("matrix index is null"); |
500
|
5122 else if (args(1).is_undefined ()) |
492
|
5123 ::error ("matrix index is a null expression"); |
|
5124 else |
500
|
5125 retval = do_matrix_index (args(1)); |
492
|
5126 break; |
|
5127 case 3: |
500
|
5128 if (args.length () <= 0) |
492
|
5129 ::error ("matrix indices are null"); |
500
|
5130 else if (args(1).is_undefined ()) |
492
|
5131 ::error ("first matrix index is a null expression"); |
500
|
5132 else if (args(2).is_undefined ()) |
492
|
5133 ::error ("second matrix index is a null expression"); |
|
5134 else |
500
|
5135 retval = do_matrix_index (args(1), args(2)); |
492
|
5136 break; |
|
5137 default: |
|
5138 ::error ("too many indices for matrix expression"); |
|
5139 break; |
|
5140 } |
|
5141 |
|
5142 return retval; |
|
5143 } |
|
5144 |
|
5145 tree_constant |
|
5146 tree_constant_rep::do_matrix_index (const tree_constant& i_arg) const |
|
5147 { |
|
5148 tree_constant retval; |
|
5149 |
|
5150 int nr = rows (); |
|
5151 int nc = columns (); |
|
5152 |
|
5153 if (user_pref.do_fortran_indexing) |
|
5154 retval = fortran_style_matrix_index (i_arg); |
|
5155 else if (nr <= 1 || nc <= 1) |
|
5156 retval = do_vector_index (i_arg); |
|
5157 else |
|
5158 ::error ("single index only valid for row or column vector"); |
|
5159 |
|
5160 return retval; |
|
5161 } |
|
5162 |
|
5163 tree_constant |
|
5164 tree_constant_rep::fortran_style_matrix_index |
|
5165 (const tree_constant& i_arg) const |
|
5166 { |
|
5167 tree_constant retval; |
|
5168 |
|
5169 tree_constant tmp_i = i_arg.make_numeric_or_magic (); |
|
5170 |
|
5171 tree_constant_rep::constant_type itype = tmp_i.const_type (); |
|
5172 |
|
5173 int nr = rows (); |
|
5174 int nc = columns (); |
|
5175 |
|
5176 switch (itype) |
|
5177 { |
|
5178 case complex_scalar_constant: |
|
5179 case scalar_constant: |
|
5180 { |
|
5181 int i = NINT (tmp_i.double_value ()); |
|
5182 int ii = fortran_row (i, nr) - 1; |
|
5183 int jj = fortran_column (i, nr) - 1; |
|
5184 if (index_check (i-1, "") < 0) |
|
5185 return tree_constant (); |
|
5186 if (range_max_check (i-1, nr * nc) < 0) |
|
5187 return tree_constant (); |
|
5188 retval = do_matrix_index (ii, jj); |
|
5189 } |
|
5190 break; |
|
5191 case complex_matrix_constant: |
|
5192 case matrix_constant: |
|
5193 { |
|
5194 Matrix mi = tmp_i.matrix_value (); |
|
5195 if (mi.rows () == 0 || mi.columns () == 0) |
|
5196 { |
|
5197 Matrix mtmp; |
|
5198 retval = tree_constant (mtmp); |
|
5199 } |
|
5200 else |
|
5201 { |
|
5202 // Yes, we really do want to call this with mi. |
|
5203 retval = fortran_style_matrix_index (mi); |
|
5204 } |
|
5205 } |
|
5206 break; |
|
5207 case string_constant: |
|
5208 gripe_string_invalid (); |
|
5209 break; |
|
5210 case range_constant: |
|
5211 gripe_range_invalid (); |
|
5212 break; |
|
5213 case magic_colon: |
|
5214 retval = do_matrix_index (magic_colon); |
|
5215 break; |
|
5216 default: |
|
5217 panic_impossible (); |
|
5218 break; |
|
5219 } |
|
5220 |
|
5221 return retval; |
|
5222 } |
|
5223 |
|
5224 tree_constant |
|
5225 tree_constant_rep::fortran_style_matrix_index (const Matrix& mi) const |
|
5226 { |
|
5227 assert (is_matrix_type ()); |
|
5228 |
|
5229 tree_constant retval; |
|
5230 |
|
5231 int nr = rows (); |
|
5232 int nc = columns (); |
|
5233 |
|
5234 int len = nr * nc; |
|
5235 |
|
5236 int index_nr = mi.rows (); |
|
5237 int index_nc = mi.columns (); |
|
5238 |
|
5239 if (index_nr >= 1 && index_nc >= 1) |
|
5240 { |
|
5241 const double *cop_out = (const double *) NULL; |
|
5242 const Complex *c_cop_out = (const Complex *) NULL; |
|
5243 int real_type = type_tag == matrix_constant; |
|
5244 if (real_type) |
|
5245 cop_out = matrix->data (); |
|
5246 else |
|
5247 c_cop_out = complex_matrix->data (); |
|
5248 |
|
5249 const double *cop_out_index = mi.data (); |
|
5250 |
|
5251 idx_vector iv (mi, 1, "", len); |
|
5252 if (! iv) |
|
5253 return tree_constant (); |
|
5254 |
|
5255 int result_size = iv.length (); |
|
5256 |
|
5257 if (nc == 1 || (nr != 1 && iv.one_zero_only ())) |
|
5258 { |
|
5259 CRMATRIX (m, cm, result_size, 1); |
|
5260 |
|
5261 for (int i = 0; i < result_size; i++) |
|
5262 { |
|
5263 int idx = iv.elem (i); |
|
5264 CRMATRIX_ASSIGN_ELEM (m, cm, i, 0, cop_out [idx], |
|
5265 c_cop_out [idx], real_type); |
|
5266 } |
|
5267 |
|
5268 ASSIGN_CRMATRIX_TO (retval, m, cm); |
|
5269 } |
|
5270 else if (nr == 1) |
|
5271 { |
|
5272 CRMATRIX (m, cm, 1, result_size); |
|
5273 |
|
5274 for (int i = 0; i < result_size; i++) |
|
5275 { |
|
5276 int idx = iv.elem (i); |
|
5277 CRMATRIX_ASSIGN_ELEM (m, cm, 0, i, cop_out [idx], |
|
5278 c_cop_out [idx], real_type); |
|
5279 } |
|
5280 |
|
5281 ASSIGN_CRMATRIX_TO (retval, m, cm); |
|
5282 } |
|
5283 else |
|
5284 { |
|
5285 CRMATRIX (m, cm, index_nr, index_nc); |
|
5286 |
|
5287 for (int j = 0; j < index_nc; j++) |
|
5288 for (int i = 0; i < index_nr; i++) |
|
5289 { |
|
5290 double tmp = *cop_out_index++; |
|
5291 int idx = tree_to_mat_idx (tmp); |
|
5292 CRMATRIX_ASSIGN_ELEM (m, cm, i, j, cop_out [idx], |
|
5293 c_cop_out [idx], real_type); |
|
5294 } |
|
5295 |
|
5296 ASSIGN_CRMATRIX_TO (retval, m, cm); |
|
5297 } |
|
5298 } |
|
5299 else |
|
5300 { |
|
5301 if (index_nr == 0 || index_nc == 0) |
|
5302 ::error ("empty matrix invalid as index"); |
|
5303 else |
|
5304 ::error ("invalid matrix index"); |
|
5305 return tree_constant (); |
|
5306 } |
|
5307 |
|
5308 return retval; |
|
5309 } |
|
5310 |
|
5311 tree_constant |
|
5312 tree_constant_rep::do_vector_index (const tree_constant& i_arg) const |
|
5313 { |
|
5314 tree_constant retval; |
|
5315 |
|
5316 tree_constant tmp_i = i_arg.make_numeric_or_range_or_magic (); |
|
5317 |
|
5318 tree_constant_rep::constant_type itype = tmp_i.const_type (); |
|
5319 |
|
5320 int nr = rows (); |
|
5321 int nc = columns (); |
|
5322 |
|
5323 int len = MAX (nr, nc); |
|
5324 |
|
5325 assert ((nr == 1 || nc == 1) && ! user_pref.do_fortran_indexing); |
|
5326 |
|
5327 int swap_indices = (nr == 1); |
|
5328 |
|
5329 switch (itype) |
|
5330 { |
|
5331 case complex_scalar_constant: |
|
5332 case scalar_constant: |
|
5333 { |
|
5334 int i = tree_to_mat_idx (tmp_i.double_value ()); |
|
5335 if (index_check (i, "") < 0) |
|
5336 return tree_constant (); |
|
5337 if (swap_indices) |
|
5338 { |
|
5339 if (range_max_check (i, nc) < 0) |
|
5340 return tree_constant (); |
|
5341 retval = do_matrix_index (0, i); |
|
5342 } |
|
5343 else |
|
5344 { |
|
5345 if (range_max_check (i, nr) < 0) |
|
5346 return tree_constant (); |
|
5347 retval = do_matrix_index (i, 0); |
|
5348 } |
|
5349 } |
|
5350 break; |
|
5351 case complex_matrix_constant: |
|
5352 case matrix_constant: |
|
5353 { |
|
5354 Matrix mi = tmp_i.matrix_value (); |
|
5355 if (mi.rows () == 0 || mi.columns () == 0) |
|
5356 { |
|
5357 Matrix mtmp; |
|
5358 retval = tree_constant (mtmp); |
|
5359 } |
|
5360 else |
|
5361 { |
|
5362 idx_vector iv (mi, user_pref.do_fortran_indexing, "", len); |
|
5363 if (! iv) |
|
5364 return tree_constant (); |
|
5365 |
|
5366 if (swap_indices) |
|
5367 { |
|
5368 if (range_max_check (iv.max (), nc) < 0) |
|
5369 return tree_constant (); |
|
5370 retval = do_matrix_index (0, iv); |
|
5371 } |
|
5372 else |
|
5373 { |
|
5374 if (range_max_check (iv.max (), nr) < 0) |
|
5375 return tree_constant (); |
|
5376 retval = do_matrix_index (iv, 0); |
|
5377 } |
|
5378 } |
|
5379 } |
|
5380 break; |
|
5381 case string_constant: |
|
5382 gripe_string_invalid (); |
|
5383 break; |
|
5384 case range_constant: |
|
5385 { |
|
5386 Range ri = tmp_i.range_value (); |
|
5387 if (len == 2 && is_zero_one (ri)) |
|
5388 { |
|
5389 if (swap_indices) |
|
5390 retval = do_matrix_index (0, 1); |
|
5391 else |
|
5392 retval = do_matrix_index (1, 0); |
|
5393 } |
|
5394 else if (len == 2 && is_one_zero (ri)) |
|
5395 { |
|
5396 retval = do_matrix_index (0, 0); |
|
5397 } |
|
5398 else |
|
5399 { |
|
5400 if (index_check (ri, "") < 0) |
|
5401 return tree_constant (); |
|
5402 if (swap_indices) |
|
5403 { |
|
5404 if (range_max_check (tree_to_mat_idx (ri.max ()), nc) < 0) |
|
5405 return tree_constant (); |
|
5406 retval = do_matrix_index (0, ri); |
|
5407 } |
|
5408 else |
|
5409 { |
|
5410 if (range_max_check (tree_to_mat_idx (ri.max ()), nr) < 0) |
|
5411 return tree_constant (); |
|
5412 retval = do_matrix_index (ri, 0); |
|
5413 } |
|
5414 } |
|
5415 } |
|
5416 break; |
|
5417 case magic_colon: |
|
5418 if (swap_indices) |
|
5419 retval = do_matrix_index (0, magic_colon); |
|
5420 else |
|
5421 retval = do_matrix_index (magic_colon, 0); |
|
5422 break; |
|
5423 default: |
|
5424 panic_impossible (); |
|
5425 break; |
|
5426 } |
|
5427 |
|
5428 return retval; |
|
5429 } |
|
5430 |
|
5431 tree_constant |
|
5432 tree_constant_rep::do_matrix_index (const tree_constant& i_arg, |
|
5433 const tree_constant& j_arg) const |
|
5434 { |
|
5435 tree_constant retval; |
|
5436 |
|
5437 tree_constant tmp_i = i_arg.make_numeric_or_range_or_magic (); |
|
5438 |
|
5439 tree_constant_rep::constant_type itype = tmp_i.const_type (); |
|
5440 |
|
5441 switch (itype) |
|
5442 { |
|
5443 case complex_scalar_constant: |
|
5444 case scalar_constant: |
|
5445 { |
|
5446 int i = tree_to_mat_idx (tmp_i.double_value ()); |
|
5447 if (index_check (i, "row") < 0) |
|
5448 return tree_constant (); |
|
5449 retval = do_matrix_index (i, j_arg); |
|
5450 } |
|
5451 break; |
|
5452 case complex_matrix_constant: |
|
5453 case matrix_constant: |
|
5454 { |
|
5455 Matrix mi = tmp_i.matrix_value (); |
|
5456 idx_vector iv (mi, user_pref.do_fortran_indexing, "row", rows ()); |
|
5457 if (! iv) |
|
5458 return tree_constant (); |
|
5459 |
|
5460 if (iv.length () == 0) |
|
5461 { |
|
5462 Matrix mtmp; |
|
5463 retval = tree_constant (mtmp); |
|
5464 } |
|
5465 else |
|
5466 retval = do_matrix_index (iv, j_arg); |
|
5467 } |
|
5468 break; |
|
5469 case string_constant: |
|
5470 gripe_string_invalid (); |
|
5471 break; |
|
5472 case range_constant: |
|
5473 { |
|
5474 Range ri = tmp_i.range_value (); |
|
5475 int nr = rows (); |
|
5476 if (nr == 2 && is_zero_one (ri)) |
|
5477 { |
|
5478 retval = do_matrix_index (1, j_arg); |
|
5479 } |
|
5480 else if (nr == 2 && is_one_zero (ri)) |
|
5481 { |
|
5482 retval = do_matrix_index (0, j_arg); |
|
5483 } |
|
5484 else |
|
5485 { |
|
5486 if (index_check (ri, "row") < 0) |
|
5487 return tree_constant (); |
|
5488 retval = do_matrix_index (ri, j_arg); |
|
5489 } |
|
5490 } |
|
5491 break; |
|
5492 case magic_colon: |
|
5493 retval = do_matrix_index (magic_colon, j_arg); |
|
5494 break; |
|
5495 default: |
|
5496 panic_impossible (); |
|
5497 break; |
|
5498 } |
|
5499 |
|
5500 return retval; |
|
5501 } |
|
5502 |
|
5503 tree_constant |
|
5504 tree_constant_rep::do_matrix_index (int i, const tree_constant& j_arg) const |
|
5505 { |
|
5506 tree_constant retval; |
|
5507 |
|
5508 tree_constant tmp_j = j_arg.make_numeric_or_range_or_magic (); |
|
5509 |
|
5510 tree_constant_rep::constant_type jtype = tmp_j.const_type (); |
|
5511 |
|
5512 int nr = rows (); |
|
5513 int nc = columns (); |
|
5514 |
|
5515 switch (jtype) |
|
5516 { |
|
5517 case complex_scalar_constant: |
|
5518 case scalar_constant: |
|
5519 { |
|
5520 int j = tree_to_mat_idx (tmp_j.double_value ()); |
|
5521 if (index_check (j, "column") < 0) |
|
5522 return tree_constant (); |
|
5523 if (range_max_check (i, j, nr, nc) < 0) |
|
5524 return tree_constant (); |
|
5525 retval = do_matrix_index (i, j); |
|
5526 } |
|
5527 break; |
|
5528 case complex_matrix_constant: |
|
5529 case matrix_constant: |
|
5530 { |
|
5531 Matrix mj = tmp_j.matrix_value (); |
|
5532 idx_vector jv (mj, user_pref.do_fortran_indexing, "column", nc); |
|
5533 if (! jv) |
|
5534 return tree_constant (); |
|
5535 |
|
5536 if (jv.length () == 0) |
|
5537 { |
|
5538 Matrix mtmp; |
|
5539 retval = tree_constant (mtmp); |
|
5540 } |
|
5541 else |
|
5542 { |
|
5543 if (range_max_check (i, jv.max (), nr, nc) < 0) |
|
5544 return tree_constant (); |
|
5545 retval = do_matrix_index (i, jv); |
|
5546 } |
|
5547 } |
|
5548 break; |
|
5549 case string_constant: |
|
5550 gripe_string_invalid (); |
|
5551 break; |
|
5552 case range_constant: |
|
5553 { |
|
5554 Range rj = tmp_j.range_value (); |
|
5555 if (nc == 2 && is_zero_one (rj)) |
|
5556 { |
|
5557 retval = do_matrix_index (i, 1); |
|
5558 } |
|
5559 else if (nc == 2 && is_one_zero (rj)) |
|
5560 { |
|
5561 retval = do_matrix_index (i, 0); |
|
5562 } |
|
5563 else |
|
5564 { |
|
5565 if (index_check (rj, "column") < 0) |
|
5566 return tree_constant (); |
|
5567 if (range_max_check (i, tree_to_mat_idx (rj.max ()), nr, nc) < 0) |
|
5568 return tree_constant (); |
|
5569 retval = do_matrix_index (i, rj); |
|
5570 } |
|
5571 } |
|
5572 break; |
|
5573 case magic_colon: |
|
5574 if (range_max_check (i, 0, nr, nc) < 0) |
|
5575 return tree_constant (); |
|
5576 retval = do_matrix_index (i, magic_colon); |
|
5577 break; |
|
5578 default: |
|
5579 panic_impossible (); |
|
5580 break; |
|
5581 } |
|
5582 |
|
5583 return retval; |
|
5584 } |
|
5585 |
|
5586 tree_constant |
|
5587 tree_constant_rep::do_matrix_index (const idx_vector& iv, |
|
5588 const tree_constant& j_arg) const |
|
5589 { |
|
5590 tree_constant retval; |
|
5591 |
|
5592 tree_constant tmp_j = j_arg.make_numeric_or_range_or_magic (); |
|
5593 |
|
5594 tree_constant_rep::constant_type jtype = tmp_j.const_type (); |
|
5595 |
|
5596 int nr = rows (); |
|
5597 int nc = columns (); |
|
5598 |
|
5599 switch (jtype) |
|
5600 { |
|
5601 case complex_scalar_constant: |
|
5602 case scalar_constant: |
|
5603 { |
|
5604 int j = tree_to_mat_idx (tmp_j.double_value ()); |
|
5605 if (index_check (j, "column") < 0) |
|
5606 return tree_constant (); |
|
5607 if (range_max_check (iv.max (), j, nr, nc) < 0) |
|
5608 return tree_constant (); |
|
5609 retval = do_matrix_index (iv, j); |
|
5610 } |
|
5611 break; |
|
5612 case complex_matrix_constant: |
|
5613 case matrix_constant: |
|
5614 { |
|
5615 Matrix mj = tmp_j.matrix_value (); |
|
5616 idx_vector jv (mj, user_pref.do_fortran_indexing, "column", nc); |
|
5617 if (! jv) |
|
5618 return tree_constant (); |
|
5619 |
|
5620 if (jv.length () == 0) |
|
5621 { |
|
5622 Matrix mtmp; |
|
5623 retval = tree_constant (mtmp); |
|
5624 } |
|
5625 else |
|
5626 { |
|
5627 if (range_max_check (iv.max (), jv.max (), nr, nc) < 0) |
|
5628 return tree_constant (); |
|
5629 retval = do_matrix_index (iv, jv); |
|
5630 } |
|
5631 } |
|
5632 break; |
|
5633 case string_constant: |
|
5634 gripe_string_invalid (); |
|
5635 break; |
|
5636 case range_constant: |
|
5637 { |
|
5638 Range rj = tmp_j.range_value (); |
|
5639 if (nc == 2 && is_zero_one (rj)) |
|
5640 { |
|
5641 retval = do_matrix_index (iv, 1); |
|
5642 } |
|
5643 else if (nc == 2 && is_one_zero (rj)) |
|
5644 { |
|
5645 retval = do_matrix_index (iv, 0); |
|
5646 } |
|
5647 else |
|
5648 { |
|
5649 if (index_check (rj, "column") < 0) |
|
5650 return tree_constant (); |
|
5651 if (range_max_check (iv.max (), tree_to_mat_idx (rj.max ()), |
|
5652 nr, nc) < 0) |
|
5653 return tree_constant (); |
|
5654 retval = do_matrix_index (iv, rj); |
|
5655 } |
|
5656 } |
|
5657 break; |
|
5658 case magic_colon: |
|
5659 if (range_max_check (iv.max (), 0, nr, nc) < 0) |
|
5660 return tree_constant (); |
|
5661 retval = do_matrix_index (iv, magic_colon); |
|
5662 break; |
|
5663 default: |
|
5664 panic_impossible (); |
|
5665 break; |
|
5666 } |
|
5667 |
|
5668 return retval; |
|
5669 } |
|
5670 |
|
5671 tree_constant |
|
5672 tree_constant_rep::do_matrix_index (const Range& ri, |
|
5673 const tree_constant& j_arg) const |
|
5674 { |
|
5675 tree_constant retval; |
|
5676 |
|
5677 tree_constant tmp_j = j_arg.make_numeric_or_range_or_magic (); |
|
5678 |
|
5679 tree_constant_rep::constant_type jtype = tmp_j.const_type (); |
|
5680 |
|
5681 int nr = rows (); |
|
5682 int nc = columns (); |
|
5683 |
|
5684 switch (jtype) |
|
5685 { |
|
5686 case complex_scalar_constant: |
|
5687 case scalar_constant: |
|
5688 { |
|
5689 int j = tree_to_mat_idx (tmp_j.double_value ()); |
|
5690 if (index_check (j, "column") < 0) |
|
5691 return tree_constant (); |
|
5692 if (range_max_check (tree_to_mat_idx (ri.max ()), j, nr, nc) < 0) |
|
5693 return tree_constant (); |
|
5694 retval = do_matrix_index (ri, j); |
|
5695 } |
|
5696 break; |
|
5697 case complex_matrix_constant: |
|
5698 case matrix_constant: |
|
5699 { |
|
5700 Matrix mj = tmp_j.matrix_value (); |
|
5701 idx_vector jv (mj, user_pref.do_fortran_indexing, "column", nc); |
|
5702 if (! jv) |
|
5703 return tree_constant (); |
|
5704 |
|
5705 if (jv.length () == 0) |
|
5706 { |
|
5707 Matrix mtmp; |
|
5708 retval = tree_constant (mtmp); |
|
5709 } |
|
5710 else |
|
5711 { |
|
5712 if (range_max_check (tree_to_mat_idx (ri.max ()), |
|
5713 jv.max (), nr, nc) < 0) |
|
5714 return tree_constant (); |
|
5715 retval = do_matrix_index (ri, jv); |
|
5716 } |
|
5717 } |
|
5718 break; |
|
5719 case string_constant: |
|
5720 gripe_string_invalid (); |
|
5721 break; |
|
5722 case range_constant: |
|
5723 { |
|
5724 Range rj = tmp_j.range_value (); |
|
5725 if (nc == 2 && is_zero_one (rj)) |
|
5726 { |
|
5727 retval = do_matrix_index (ri, 1); |
|
5728 } |
|
5729 else if (nc == 2 && is_one_zero (rj)) |
|
5730 { |
|
5731 retval = do_matrix_index (ri, 0); |
|
5732 } |
|
5733 else |
|
5734 { |
|
5735 if (index_check (rj, "column") < 0) |
|
5736 return tree_constant (); |
|
5737 if (range_max_check (tree_to_mat_idx (ri.max ()), |
|
5738 tree_to_mat_idx (rj.max ()), nr, nc) < 0) |
|
5739 return tree_constant (); |
|
5740 retval = do_matrix_index (ri, rj); |
|
5741 } |
|
5742 } |
|
5743 break; |
|
5744 case magic_colon: |
|
5745 retval = do_matrix_index (ri, magic_colon); |
|
5746 break; |
|
5747 default: |
|
5748 panic_impossible (); |
|
5749 break; |
|
5750 } |
|
5751 |
|
5752 return retval; |
|
5753 } |
|
5754 |
|
5755 tree_constant |
|
5756 tree_constant_rep::do_matrix_index (tree_constant_rep::constant_type mci, |
|
5757 const tree_constant& j_arg) const |
|
5758 { |
|
5759 tree_constant retval; |
|
5760 |
|
5761 tree_constant tmp_j = j_arg.make_numeric_or_range_or_magic (); |
|
5762 |
|
5763 tree_constant_rep::constant_type jtype = tmp_j.const_type (); |
|
5764 |
|
5765 int nr = rows (); |
|
5766 int nc = columns (); |
|
5767 |
|
5768 switch (jtype) |
|
5769 { |
|
5770 case complex_scalar_constant: |
|
5771 case scalar_constant: |
|
5772 { |
|
5773 int j = tree_to_mat_idx (tmp_j.double_value ()); |
|
5774 if (index_check (j, "column") < 0) |
|
5775 return tree_constant (); |
|
5776 if (range_max_check (0, j, nr, nc) < 0) |
|
5777 return tree_constant (); |
|
5778 retval = do_matrix_index (magic_colon, j); |
|
5779 } |
|
5780 break; |
|
5781 case complex_matrix_constant: |
|
5782 case matrix_constant: |
|
5783 { |
|
5784 Matrix mj = tmp_j.matrix_value (); |
|
5785 idx_vector jv (mj, user_pref.do_fortran_indexing, "column", nc); |
|
5786 if (! jv) |
|
5787 return tree_constant (); |
|
5788 |
|
5789 if (jv.length () == 0) |
|
5790 { |
|
5791 Matrix mtmp; |
|
5792 retval = tree_constant (mtmp); |
|
5793 } |
|
5794 else |
|
5795 { |
|
5796 if (range_max_check (0, jv.max (), nr, nc) < 0) |
|
5797 return tree_constant (); |
|
5798 retval = do_matrix_index (magic_colon, jv); |
|
5799 } |
|
5800 } |
|
5801 break; |
|
5802 case string_constant: |
|
5803 gripe_string_invalid (); |
|
5804 break; |
|
5805 case range_constant: |
|
5806 { |
|
5807 Range rj = tmp_j.range_value (); |
|
5808 if (nc == 2 && is_zero_one (rj)) |
|
5809 { |
|
5810 retval = do_matrix_index (magic_colon, 1); |
|
5811 } |
|
5812 else if (nc == 2 && is_one_zero (rj)) |
|
5813 { |
|
5814 retval = do_matrix_index (magic_colon, 0); |
|
5815 } |
|
5816 else |
|
5817 { |
|
5818 if (index_check (rj, "column") < 0) |
|
5819 return tree_constant (); |
|
5820 if (range_max_check (0, tree_to_mat_idx (rj.max ()), nr, nc) < 0) |
|
5821 return tree_constant (); |
|
5822 retval = do_matrix_index (magic_colon, rj); |
|
5823 } |
|
5824 } |
|
5825 break; |
|
5826 case magic_colon: |
|
5827 retval = do_matrix_index (magic_colon, magic_colon); |
|
5828 break; |
|
5829 default: |
|
5830 panic_impossible (); |
|
5831 break; |
|
5832 } |
|
5833 |
|
5834 return retval; |
|
5835 } |
|
5836 |
|
5837 tree_constant |
|
5838 tree_constant_rep::do_matrix_index (int i, int j) const |
|
5839 { |
|
5840 tree_constant retval; |
|
5841 |
|
5842 if (type_tag == matrix_constant) |
|
5843 retval = tree_constant (matrix->elem (i, j)); |
|
5844 else |
|
5845 retval = tree_constant (complex_matrix->elem (i, j)); |
|
5846 |
|
5847 return retval; |
|
5848 } |
|
5849 |
|
5850 tree_constant |
|
5851 tree_constant_rep::do_matrix_index (int i, const idx_vector& jv) const |
|
5852 { |
|
5853 tree_constant retval; |
|
5854 |
|
5855 int jlen = jv.capacity (); |
|
5856 |
|
5857 CRMATRIX (m, cm, 1, jlen); |
|
5858 |
|
5859 for (int j = 0; j < jlen; j++) |
|
5860 { |
|
5861 int col = jv.elem (j); |
|
5862 CRMATRIX_ASSIGN_REP_ELEM (m, cm, 0, j, i, col); |
|
5863 } |
|
5864 ASSIGN_CRMATRIX_TO (retval, m, cm); |
|
5865 |
|
5866 return retval; |
|
5867 } |
|
5868 |
|
5869 tree_constant |
|
5870 tree_constant_rep::do_matrix_index (int i, const Range& rj) const |
|
5871 { |
|
5872 tree_constant retval; |
|
5873 |
|
5874 int jlen = rj.nelem (); |
|
5875 |
|
5876 CRMATRIX (m, cm, 1, jlen); |
|
5877 |
|
5878 double b = rj.base (); |
|
5879 double increment = rj.inc (); |
|
5880 for (int j = 0; j < jlen; j++) |
|
5881 { |
|
5882 double tmp = b + j * increment; |
|
5883 int col = tree_to_mat_idx (tmp); |
|
5884 CRMATRIX_ASSIGN_REP_ELEM (m, cm, 0, j, i, col); |
|
5885 } |
|
5886 |
|
5887 ASSIGN_CRMATRIX_TO (retval, m, cm); |
|
5888 |
|
5889 return retval; |
|
5890 } |
|
5891 |
|
5892 tree_constant |
|
5893 tree_constant_rep::do_matrix_index |
|
5894 (int i, tree_constant_rep::constant_type mcj) const |
|
5895 { |
|
5896 assert (mcj == magic_colon); |
|
5897 |
|
5898 tree_constant retval; |
|
5899 |
|
5900 int nc = columns (); |
|
5901 |
|
5902 CRMATRIX (m, cm, 1, nc); |
|
5903 |
|
5904 for (int j = 0; j < nc; j++) |
|
5905 { |
|
5906 CRMATRIX_ASSIGN_REP_ELEM (m, cm, 0, j, i, j); |
|
5907 } |
|
5908 |
|
5909 ASSIGN_CRMATRIX_TO (retval, m, cm); |
|
5910 |
|
5911 return retval; |
|
5912 } |
|
5913 |
|
5914 tree_constant |
|
5915 tree_constant_rep::do_matrix_index (const idx_vector& iv, int j) const |
|
5916 { |
|
5917 tree_constant retval; |
|
5918 |
|
5919 int ilen = iv.capacity (); |
|
5920 |
|
5921 CRMATRIX (m, cm, ilen, 1); |
|
5922 |
|
5923 for (int i = 0; i < ilen; i++) |
|
5924 { |
|
5925 int row = iv.elem (i); |
|
5926 CRMATRIX_ASSIGN_REP_ELEM (m, cm, i, 0, row, j); |
|
5927 } |
|
5928 |
|
5929 ASSIGN_CRMATRIX_TO (retval, m, cm); |
|
5930 |
|
5931 return retval; |
|
5932 } |
|
5933 |
|
5934 tree_constant |
|
5935 tree_constant_rep::do_matrix_index (const idx_vector& iv, |
|
5936 const idx_vector& jv) const |
|
5937 { |
|
5938 tree_constant retval; |
|
5939 |
|
5940 int ilen = iv.capacity (); |
|
5941 int jlen = jv.capacity (); |
|
5942 |
|
5943 CRMATRIX (m, cm, ilen, jlen); |
|
5944 |
|
5945 for (int i = 0; i < ilen; i++) |
|
5946 { |
|
5947 int row = iv.elem (i); |
|
5948 for (int j = 0; j < jlen; j++) |
|
5949 { |
|
5950 int col = jv.elem (j); |
|
5951 CRMATRIX_ASSIGN_REP_ELEM (m, cm, i, j, row, col); |
|
5952 } |
|
5953 } |
|
5954 |
|
5955 ASSIGN_CRMATRIX_TO (retval, m, cm); |
|
5956 |
|
5957 return retval; |
|
5958 } |
|
5959 |
|
5960 tree_constant |
|
5961 tree_constant_rep::do_matrix_index (const idx_vector& iv, |
|
5962 const Range& rj) const |
|
5963 { |
|
5964 tree_constant retval; |
|
5965 |
|
5966 int ilen = iv.capacity (); |
|
5967 int jlen = rj.nelem (); |
|
5968 |
|
5969 CRMATRIX (m, cm, ilen, jlen); |
|
5970 |
|
5971 double b = rj.base (); |
|
5972 double increment = rj.inc (); |
|
5973 |
|
5974 for (int i = 0; i < ilen; i++) |
|
5975 { |
|
5976 int row = iv.elem (i); |
|
5977 for (int j = 0; j < jlen; j++) |
|
5978 { |
|
5979 double tmp = b + j * increment; |
|
5980 int col = tree_to_mat_idx (tmp); |
|
5981 CRMATRIX_ASSIGN_REP_ELEM (m, cm, i, j, row, col); |
|
5982 } |
|
5983 } |
|
5984 |
|
5985 ASSIGN_CRMATRIX_TO (retval, m, cm); |
|
5986 |
|
5987 return retval; |
|
5988 } |
|
5989 |
|
5990 tree_constant |
|
5991 tree_constant_rep::do_matrix_index |
|
5992 (const idx_vector& iv, tree_constant_rep::constant_type mcj) const |
|
5993 { |
|
5994 assert (mcj == magic_colon); |
|
5995 |
|
5996 tree_constant retval; |
|
5997 |
|
5998 int nc = columns (); |
|
5999 int ilen = iv.capacity (); |
|
6000 |
|
6001 CRMATRIX (m, cm, ilen, nc); |
|
6002 |
|
6003 for (int j = 0; j < nc; j++) |
|
6004 { |
|
6005 for (int i = 0; i < ilen; i++) |
|
6006 { |
|
6007 int row = iv.elem (i); |
|
6008 CRMATRIX_ASSIGN_REP_ELEM (m, cm, i, j, row, j); |
|
6009 } |
|
6010 } |
|
6011 |
|
6012 ASSIGN_CRMATRIX_TO (retval, m, cm); |
|
6013 |
|
6014 return retval; |
|
6015 } |
|
6016 |
|
6017 tree_constant |
|
6018 tree_constant_rep::do_matrix_index (const Range& ri, int j) const |
|
6019 { |
|
6020 tree_constant retval; |
|
6021 |
|
6022 int ilen = ri.nelem (); |
|
6023 |
|
6024 CRMATRIX (m, cm, ilen, 1); |
|
6025 |
|
6026 double b = ri.base (); |
|
6027 double increment = ri.inc (); |
|
6028 for (int i = 0; i < ilen; i++) |
|
6029 { |
|
6030 double tmp = b + i * increment; |
|
6031 int row = tree_to_mat_idx (tmp); |
|
6032 CRMATRIX_ASSIGN_REP_ELEM (m, cm, i, 0, row, j); |
|
6033 } |
|
6034 |
|
6035 ASSIGN_CRMATRIX_TO (retval, m, cm); |
|
6036 |
|
6037 return retval; |
|
6038 } |
|
6039 |
|
6040 tree_constant |
|
6041 tree_constant_rep::do_matrix_index (const Range& ri, |
|
6042 const idx_vector& jv) const |
|
6043 { |
|
6044 tree_constant retval; |
|
6045 |
|
6046 int ilen = ri.nelem (); |
|
6047 int jlen = jv.capacity (); |
|
6048 |
|
6049 CRMATRIX (m, cm, ilen, jlen); |
|
6050 |
|
6051 double b = ri.base (); |
|
6052 double increment = ri.inc (); |
|
6053 for (int i = 0; i < ilen; i++) |
|
6054 { |
|
6055 double tmp = b + i * increment; |
|
6056 int row = tree_to_mat_idx (tmp); |
|
6057 for (int j = 0; j < jlen; j++) |
|
6058 { |
|
6059 int col = jv.elem (j); |
|
6060 CRMATRIX_ASSIGN_REP_ELEM (m, cm, i, j, row, col); |
|
6061 } |
|
6062 } |
|
6063 |
|
6064 ASSIGN_CRMATRIX_TO (retval, m, cm); |
|
6065 |
|
6066 return retval; |
|
6067 } |
|
6068 |
|
6069 tree_constant |
|
6070 tree_constant_rep::do_matrix_index (const Range& ri, const Range& rj) const |
|
6071 { |
|
6072 tree_constant retval; |
|
6073 |
|
6074 int ilen = ri.nelem (); |
|
6075 int jlen = rj.nelem (); |
|
6076 |
|
6077 CRMATRIX (m, cm, ilen, jlen); |
|
6078 |
|
6079 double ib = ri.base (); |
|
6080 double iinc = ri.inc (); |
|
6081 double jb = rj.base (); |
|
6082 double jinc = rj.inc (); |
|
6083 |
|
6084 for (int i = 0; i < ilen; i++) |
|
6085 { |
|
6086 double itmp = ib + i * iinc; |
|
6087 int row = tree_to_mat_idx (itmp); |
|
6088 for (int j = 0; j < jlen; j++) |
|
6089 { |
|
6090 double jtmp = jb + j * jinc; |
|
6091 int col = tree_to_mat_idx (jtmp); |
|
6092 |
|
6093 CRMATRIX_ASSIGN_REP_ELEM (m, cm, i, j, row, col); |
|
6094 } |
|
6095 } |
|
6096 |
|
6097 ASSIGN_CRMATRIX_TO (retval, m, cm); |
|
6098 |
|
6099 return retval; |
|
6100 } |
|
6101 |
|
6102 tree_constant |
|
6103 tree_constant_rep::do_matrix_index |
|
6104 (const Range& ri, tree_constant_rep::constant_type mcj) const |
|
6105 { |
|
6106 assert (mcj == magic_colon); |
|
6107 |
|
6108 tree_constant retval; |
|
6109 |
|
6110 int nc = columns (); |
|
6111 |
|
6112 int ilen = ri.nelem (); |
|
6113 |
|
6114 CRMATRIX (m, cm, ilen, nc); |
|
6115 |
|
6116 double ib = ri.base (); |
|
6117 double iinc = ri.inc (); |
|
6118 |
|
6119 for (int i = 0; i < ilen; i++) |
|
6120 { |
|
6121 double itmp = ib + i * iinc; |
|
6122 int row = tree_to_mat_idx (itmp); |
|
6123 for (int j = 0; j < nc; j++) |
|
6124 { |
|
6125 CRMATRIX_ASSIGN_REP_ELEM (m, cm, i, j, row, j); |
|
6126 } |
|
6127 } |
|
6128 |
|
6129 ASSIGN_CRMATRIX_TO (retval, m, cm); |
|
6130 |
|
6131 return retval; |
|
6132 } |
|
6133 |
|
6134 tree_constant |
|
6135 tree_constant_rep::do_matrix_index (tree_constant_rep::constant_type mci, |
|
6136 int j) const |
|
6137 { |
|
6138 assert (mci == magic_colon); |
|
6139 |
|
6140 tree_constant retval; |
|
6141 |
|
6142 int nr = rows (); |
|
6143 |
|
6144 CRMATRIX (m, cm, nr, 1); |
|
6145 |
|
6146 for (int i = 0; i < nr; i++) |
|
6147 { |
|
6148 CRMATRIX_ASSIGN_REP_ELEM (m, cm, i, 0, i, j); |
|
6149 } |
|
6150 |
|
6151 ASSIGN_CRMATRIX_TO (retval, m, cm); |
|
6152 |
|
6153 return retval; |
|
6154 } |
|
6155 |
|
6156 tree_constant |
|
6157 tree_constant_rep::do_matrix_index (tree_constant_rep::constant_type mci, |
|
6158 const idx_vector& jv) const |
|
6159 { |
|
6160 assert (mci == magic_colon); |
|
6161 |
|
6162 tree_constant retval; |
|
6163 |
|
6164 int nr = rows (); |
|
6165 int jlen = jv.capacity (); |
|
6166 |
|
6167 CRMATRIX (m, cm, nr, jlen); |
|
6168 |
|
6169 for (int i = 0; i < nr; i++) |
|
6170 { |
|
6171 for (int j = 0; j < jlen; j++) |
|
6172 { |
|
6173 int col = jv.elem (j); |
|
6174 CRMATRIX_ASSIGN_REP_ELEM (m, cm, i, j, i, col); |
|
6175 } |
|
6176 } |
|
6177 |
|
6178 ASSIGN_CRMATRIX_TO (retval, m, cm); |
|
6179 |
|
6180 return retval; |
|
6181 } |
|
6182 |
|
6183 tree_constant |
|
6184 tree_constant_rep::do_matrix_index (tree_constant_rep::constant_type mci, |
|
6185 const Range& rj) const |
|
6186 { |
|
6187 assert (mci == magic_colon); |
|
6188 |
|
6189 tree_constant retval; |
|
6190 |
|
6191 int nr = rows (); |
|
6192 int jlen = rj.nelem (); |
|
6193 |
|
6194 CRMATRIX (m, cm, nr, jlen); |
|
6195 |
|
6196 double jb = rj.base (); |
|
6197 double jinc = rj.inc (); |
|
6198 |
|
6199 for (int j = 0; j < jlen; j++) |
|
6200 { |
|
6201 double jtmp = jb + j * jinc; |
|
6202 int col = tree_to_mat_idx (jtmp); |
|
6203 for (int i = 0; i < nr; i++) |
|
6204 { |
|
6205 CRMATRIX_ASSIGN_REP_ELEM (m, cm, i, j, i, col); |
|
6206 } |
|
6207 } |
|
6208 |
|
6209 ASSIGN_CRMATRIX_TO (retval, m, cm); |
|
6210 |
|
6211 return retval; |
|
6212 } |
|
6213 |
|
6214 tree_constant |
|
6215 tree_constant_rep::do_matrix_index (tree_constant_rep::constant_type mci, |
|
6216 tree_constant_rep::constant_type mcj) const |
|
6217 { |
|
6218 assert (mci == magic_colon && mcj == magic_colon); |
|
6219 |
|
6220 return tree_constant (*this); |
|
6221 } |
|
6222 |
|
6223 tree_constant |
|
6224 tree_constant_rep::do_matrix_index |
|
6225 (tree_constant_rep::constant_type mci) const |
|
6226 { |
|
6227 assert (mci == magic_colon); |
|
6228 |
|
6229 tree_constant retval; |
|
6230 int nr = rows (); |
|
6231 int nc = columns (); |
|
6232 int size = nr * nc; |
|
6233 if (size > 0) |
|
6234 { |
|
6235 CRMATRIX (m, cm, size, 1); |
|
6236 int idx = 0; |
|
6237 for (int j = 0; j < nc; j++) |
|
6238 for (int i = 0; i < nr; i++) |
|
6239 { |
|
6240 CRMATRIX_ASSIGN_REP_ELEM (m, cm, idx, 0, i, j); |
|
6241 idx++; |
|
6242 } |
|
6243 ASSIGN_CRMATRIX_TO (retval, m, cm); |
|
6244 } |
|
6245 return retval; |
|
6246 } |
|
6247 |
|
6248 /* |
|
6249 ;;; Local Variables: *** |
|
6250 ;;; mode: C++ *** |
|
6251 ;;; page-delimiter: "^/\\*" *** |
|
6252 ;;; End: *** |
|
6253 */ |