1
|
1 // Helper functions for arithmetic operations. -*- C++ -*- |
|
2 // Used by the tree class. |
|
3 /* |
|
4 |
390
|
5 Copyright (C) 1992, 1993, 1994 John W. Eaton |
1
|
6 |
|
7 This file is part of Octave. |
|
8 |
|
9 Octave is free software; you can redistribute it and/or modify it |
|
10 under the terms of the GNU General Public License as published by the |
|
11 Free Software Foundation; either version 2, or (at your option) any |
|
12 later version. |
|
13 |
|
14 Octave is distributed in the hope that it will be useful, but WITHOUT |
|
15 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or |
|
16 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License |
|
17 for more details. |
|
18 |
|
19 You should have received a copy of the GNU General Public License |
|
20 along with Octave; see the file COPYING. If not, write to the Free |
|
21 Software Foundation, 675 Mass Ave, Cambridge, MA 02139, USA. |
|
22 |
|
23 */ |
|
24 |
240
|
25 #ifdef HAVE_CONFIG_H |
|
26 #include "config.h" |
1
|
27 #endif |
|
28 |
|
29 #include <ctype.h> |
|
30 #include <setjmp.h> |
|
31 #include <math.h> |
164
|
32 #include <Complex.h> |
1
|
33 |
|
34 #include "error.h" |
|
35 #include "gripes.h" |
|
36 #include "utils.h" |
|
37 #include "mappers.h" |
|
38 #include "user-prefs.h" |
|
39 #include "tree-const.h" |
|
40 #include "arith-ops.h" |
|
41 #include "unwind-prot.h" |
|
42 #include "xpow.h" |
|
43 #include "xdiv.h" |
|
44 |
|
45 #if defined (HAVE_ISINF) || (defined (HAVE_FINITE) && defined (HAVE_ISNAN)) |
|
46 #define DIVIDE_BY_ZERO_ERROR \ |
|
47 do \ |
|
48 { \ |
|
49 if (user_pref.warn_divide_by_zero) \ |
|
50 warning ("division by zero"); \ |
|
51 } \ |
|
52 while (0) |
|
53 #else |
|
54 #define DIVIDE_BY_ZERO_ERROR \ |
|
55 do \ |
|
56 { \ |
|
57 error ("division by zero attempted"); \ |
|
58 return tree_constant (); \ |
|
59 } \ |
|
60 while (0) |
|
61 #endif |
|
62 |
|
63 // But first, some stupid functions that don\'t deserve to be in the |
|
64 // Matrix class... |
|
65 |
|
66 enum |
|
67 Matrix_bool_op |
|
68 { |
|
69 Matrix_LT, |
|
70 Matrix_LE, |
|
71 Matrix_EQ, |
|
72 Matrix_GE, |
|
73 Matrix_GT, |
|
74 Matrix_NE, |
|
75 Matrix_AND, |
|
76 Matrix_OR, |
|
77 }; |
|
78 |
|
79 /* |
143
|
80 * Check row and column dimensions for binary matrix operations. |
|
81 */ |
|
82 static inline int |
164
|
83 m_add_conform (const Matrix& a, const Matrix& b, int warn) |
143
|
84 { |
|
85 int ar = a.rows (); |
|
86 int ac = a.columns (); |
|
87 int br = b.rows (); |
|
88 int bc = b.columns (); |
|
89 |
|
90 int ok = (ar == br && ac == bc); |
|
91 |
|
92 if (! ok && warn) |
|
93 gripe_nonconformant (ar, ac, br, bc); |
|
94 |
|
95 return ok; |
|
96 } |
|
97 |
|
98 static inline int |
164
|
99 m_add_conform (const Matrix& a, const ComplexMatrix& b, int warn) |
143
|
100 { |
|
101 int ar = a.rows (); |
|
102 int ac = a.columns (); |
|
103 int br = b.rows (); |
|
104 int bc = b.columns (); |
|
105 |
|
106 int ok = (ar == br && ac == bc); |
|
107 |
|
108 if (! ok && warn) |
|
109 gripe_nonconformant (ar, ac, br, bc); |
|
110 |
|
111 return ok; |
|
112 } |
|
113 |
|
114 static inline int |
164
|
115 m_add_conform (const ComplexMatrix& a, const Matrix& b, int warn) |
143
|
116 { |
|
117 int ar = a.rows (); |
|
118 int ac = a.columns (); |
|
119 int br = b.rows (); |
|
120 int bc = b.columns (); |
|
121 |
|
122 int ok = (ar == br && ac == bc); |
|
123 |
|
124 if (! ok && warn) |
|
125 gripe_nonconformant (ar, ac, br, bc); |
|
126 |
|
127 return ok; |
|
128 } |
|
129 |
|
130 static inline int |
164
|
131 m_add_conform (const ComplexMatrix& a, const ComplexMatrix& b, int warn) |
143
|
132 { |
|
133 int ar = a.rows (); |
|
134 int ac = a.columns (); |
|
135 int br = b.rows (); |
|
136 int bc = b.columns (); |
|
137 |
|
138 int ok = (ar == br && ac == bc); |
|
139 |
|
140 if (! ok && warn) |
|
141 gripe_nonconformant (ar, ac, br, bc); |
|
142 |
|
143 return ok; |
|
144 } |
|
145 |
|
146 static inline int |
164
|
147 m_mul_conform (const Matrix& a, const Matrix& b, int warn) |
143
|
148 { |
|
149 int ac = a.columns (); |
|
150 int br = b.rows (); |
|
151 |
|
152 int ok = (ac == br); |
|
153 |
|
154 if (! ok && warn) |
|
155 gripe_nonconformant (a.rows (), ac, br, b.columns ()); |
|
156 |
|
157 return ok; |
|
158 } |
|
159 |
|
160 static inline int |
164
|
161 m_mul_conform (const Matrix& a, const ComplexMatrix& b, int warn) |
143
|
162 { |
|
163 int ac = a.columns (); |
|
164 int br = b.rows (); |
|
165 |
|
166 int ok = (ac == br); |
|
167 |
|
168 if (! ok && warn) |
|
169 gripe_nonconformant (a.rows (), ac, br, b.columns ()); |
|
170 |
|
171 return ok; |
|
172 } |
|
173 |
|
174 static inline int |
164
|
175 m_mul_conform (const ComplexMatrix& a, const Matrix& b, int warn) |
143
|
176 { |
|
177 int ac = a.columns (); |
|
178 int br = b.rows (); |
|
179 |
|
180 int ok = (ac == br); |
|
181 |
|
182 if (! ok && warn) |
|
183 gripe_nonconformant (a.rows (), ac, br, b.columns ()); |
|
184 |
|
185 return ok; |
|
186 } |
|
187 |
|
188 static inline int |
164
|
189 m_mul_conform (const ComplexMatrix& a, const ComplexMatrix& b, int warn) |
143
|
190 { |
|
191 int ac = a.columns (); |
|
192 int br = b.rows (); |
|
193 |
|
194 int ok = (a.columns () == br); |
|
195 |
|
196 if (! ok && warn) |
|
197 gripe_nonconformant (a.rows (), ac, br, b.columns ()); |
|
198 |
|
199 return ok; |
|
200 } |
|
201 |
|
202 /* |
1
|
203 * Stupid binary comparison operations like the ones Matlab provides. |
|
204 * One for each type combination, in the order given here: |
|
205 * |
|
206 * op2 \ op1: s m cs cm |
|
207 * +-- +---+---+----+----+ |
|
208 * scalar | | * | 3 | * | 9 | |
|
209 * +---+---+----+----+ |
|
210 * matrix | 1 | 4 | 7 | 10 | |
|
211 * +---+---+----+----+ |
|
212 * complex_scalar | * | 5 | * | 11 | |
|
213 * +---+---+----+----+ |
|
214 * complex_matrix | 2 | 6 | 8 | 12 | |
|
215 * +---+---+----+----+ |
|
216 */ |
|
217 |
|
218 /* 1 */ |
|
219 static Matrix |
164
|
220 mx_stupid_bool_op (Matrix_bool_op op, double s, const Matrix& a) |
1
|
221 { |
|
222 int ar = a.rows (); |
|
223 int ac = a.columns (); |
143
|
224 |
1
|
225 Matrix t (ar, ac); |
143
|
226 |
1
|
227 for (int j = 0; j < ac; j++) |
|
228 for (int i = 0; i < ar; i++) |
|
229 { |
|
230 switch (op) |
|
231 { |
|
232 case Matrix_LT: |
|
233 t.elem (i,j) = s < a.elem (i,j); |
|
234 break; |
|
235 case Matrix_LE: |
|
236 t.elem (i,j) = s <= a.elem (i,j); |
|
237 break; |
|
238 case Matrix_EQ: |
|
239 t.elem (i,j) = s == a.elem (i,j); |
|
240 break; |
|
241 case Matrix_GE: |
|
242 t.elem (i,j) = s >= a.elem (i,j); |
|
243 break; |
|
244 case Matrix_GT: |
|
245 t.elem (i,j) = s > a.elem (i,j); |
|
246 break; |
|
247 case Matrix_NE: |
|
248 t.elem (i,j) = s != a.elem (i,j); |
|
249 break; |
|
250 case Matrix_AND: |
|
251 t.elem (i,j) = s && a.elem (i,j); |
|
252 break; |
|
253 case Matrix_OR: |
|
254 t.elem (i,j) = s || a.elem (i,j); |
|
255 break; |
|
256 default: |
|
257 panic_impossible (); |
|
258 break; |
|
259 } |
|
260 } |
143
|
261 |
1
|
262 return t; |
|
263 } |
|
264 |
|
265 /* 2 */ |
|
266 static Matrix |
164
|
267 mx_stupid_bool_op (Matrix_bool_op op, double s, const ComplexMatrix& a) |
1
|
268 { |
|
269 int ar = a.rows (); |
|
270 int ac = a.columns (); |
143
|
271 |
1
|
272 Matrix t (ar, ac); |
143
|
273 |
1
|
274 for (int j = 0; j < ac; j++) |
|
275 for (int i = 0; i < ar; i++) |
|
276 { |
|
277 switch (op) |
|
278 { |
|
279 case Matrix_LT: |
|
280 t.elem (i,j) = s < real (a.elem (i,j)); |
|
281 break; |
|
282 case Matrix_LE: |
|
283 t.elem (i,j) = s <= real (a.elem (i,j)); |
|
284 break; |
|
285 case Matrix_EQ: |
|
286 t.elem (i,j) = s == a.elem (i,j); |
|
287 break; |
|
288 case Matrix_GE: |
|
289 t.elem (i,j) = s >= real (a.elem (i,j)); |
|
290 break; |
|
291 case Matrix_GT: |
|
292 t.elem (i,j) = s > real (a.elem (i,j)); |
|
293 break; |
|
294 case Matrix_NE: |
|
295 t.elem (i,j) = s != a.elem (i,j); |
|
296 break; |
|
297 case Matrix_AND: |
|
298 t.elem (i,j) = s && (a.elem (i,j) != 0.0); |
|
299 break; |
|
300 case Matrix_OR: |
|
301 t.elem (i,j) = s || (a.elem (i,j) != 0.0); |
|
302 break; |
|
303 default: |
|
304 panic_impossible (); |
|
305 break; |
|
306 } |
|
307 } |
143
|
308 |
1
|
309 return t; |
|
310 } |
|
311 |
|
312 /* 3 */ |
|
313 static Matrix |
164
|
314 mx_stupid_bool_op (Matrix_bool_op op, const Matrix& a, double s) |
1
|
315 { |
|
316 int ar = a.rows (); |
|
317 int ac = a.columns (); |
143
|
318 |
1
|
319 Matrix t (ar, ac); |
143
|
320 |
1
|
321 for (int j = 0; j < ac; j++) |
|
322 for (int i = 0; i < ar; i++) |
|
323 { |
|
324 switch (op) |
|
325 { |
|
326 case Matrix_LT: |
|
327 t.elem (i,j) = a.elem (i,j) < s; |
|
328 break; |
|
329 case Matrix_LE: |
|
330 t.elem (i,j) = a.elem (i,j) <= s; |
|
331 break; |
|
332 case Matrix_EQ: |
|
333 t.elem (i,j) = a.elem (i,j) == s; |
|
334 break; |
|
335 case Matrix_GE: |
|
336 t.elem (i,j) = a.elem (i,j) >= s; |
|
337 break; |
|
338 case Matrix_GT: |
|
339 t.elem (i,j) = a.elem (i,j) > s; |
|
340 break; |
|
341 case Matrix_NE: |
|
342 t.elem (i,j) = a.elem (i,j) != s; |
|
343 break; |
|
344 case Matrix_AND: |
|
345 t.elem (i,j) = a.elem (i,j) && s; |
|
346 break; |
|
347 case Matrix_OR: |
|
348 t.elem (i,j) = a.elem (i,j) || s; |
|
349 break; |
|
350 default: |
|
351 panic_impossible (); |
|
352 break; |
|
353 } |
|
354 } |
143
|
355 |
1
|
356 return t; |
|
357 } |
|
358 |
|
359 /* 4 */ |
|
360 static Matrix |
164
|
361 mx_stupid_bool_op (Matrix_bool_op op, const Matrix& a, const Complex& s) |
1
|
362 { |
|
363 int ar = a.rows (); |
|
364 int ac = a.columns (); |
143
|
365 |
1
|
366 Matrix t (ar, ac); |
143
|
367 |
1
|
368 for (int j = 0; j < ac; j++) |
|
369 for (int i = 0; i < ar; i++) |
|
370 { |
|
371 switch (op) |
|
372 { |
|
373 case Matrix_LT: |
|
374 t.elem (i,j) = a.elem (i,j) < real (s); |
|
375 break; |
|
376 case Matrix_LE: |
|
377 t.elem (i,j) = a.elem (i,j) <= real (s); |
|
378 break; |
|
379 case Matrix_EQ: |
|
380 t.elem (i,j) = a.elem (i,j) == s; |
|
381 break; |
|
382 case Matrix_GE: |
|
383 t.elem (i,j) = a.elem (i,j) >= real (s); |
|
384 break; |
|
385 case Matrix_GT: |
|
386 t.elem (i,j) = a.elem (i,j) > real (s); |
|
387 break; |
|
388 case Matrix_NE: |
|
389 t.elem (i,j) = a.elem (i,j) != s; |
|
390 break; |
|
391 case Matrix_AND: |
|
392 t.elem (i,j) = a.elem (i,j) && (s != 0.0); |
|
393 break; |
|
394 case Matrix_OR: |
|
395 t.elem (i,j) = a.elem (i,j) || (s != 0.0); |
|
396 break; |
|
397 default: |
|
398 panic_impossible (); |
|
399 break; |
|
400 } |
|
401 } |
143
|
402 |
1
|
403 return t; |
|
404 } |
|
405 |
|
406 /* 5 */ |
|
407 static Matrix |
164
|
408 mx_stupid_bool_op (Matrix_bool_op op, const Matrix& a, const Matrix& b) |
1
|
409 { |
143
|
410 if (! m_add_conform (a, b, 1)) |
|
411 return Matrix (); |
|
412 |
1
|
413 int ar = a.rows (); |
|
414 int ac = a.columns (); |
|
415 |
|
416 Matrix c (ar, ac); |
|
417 |
|
418 for (int j = 0; j < ac; j++) |
|
419 for (int i = 0; i < ar; i++) |
|
420 { |
|
421 switch (op) |
|
422 { |
|
423 case Matrix_LT: |
|
424 c.elem (i, j) = a.elem (i, j) < b.elem (i, j); |
|
425 break; |
|
426 case Matrix_LE: |
|
427 c.elem (i, j) = a.elem (i, j) <= b.elem (i, j); |
|
428 break; |
|
429 case Matrix_EQ: |
|
430 c.elem (i, j) = a.elem (i, j) == b.elem (i, j); |
|
431 break; |
|
432 case Matrix_GE: |
|
433 c.elem (i, j) = a.elem (i, j) >= b.elem (i, j); |
|
434 break; |
|
435 case Matrix_GT: |
|
436 c.elem (i, j) = a.elem (i, j) > b.elem (i, j); |
|
437 break; |
|
438 case Matrix_NE: |
|
439 c.elem (i, j) = a.elem (i, j) != b.elem (i, j); |
|
440 break; |
|
441 case Matrix_AND: |
|
442 c.elem (i, j) = a.elem (i, j) && b.elem (i, j); |
|
443 break; |
|
444 case Matrix_OR: |
|
445 c.elem (i, j) = a.elem (i, j) || b.elem (i, j); |
|
446 break; |
|
447 default: |
|
448 panic_impossible (); |
|
449 break; |
|
450 } |
|
451 } |
143
|
452 |
1
|
453 return c; |
|
454 } |
|
455 |
|
456 /* 6 */ |
|
457 static Matrix |
164
|
458 mx_stupid_bool_op (Matrix_bool_op op, const Matrix& a, const ComplexMatrix& b) |
1
|
459 { |
143
|
460 if (! m_add_conform (a, b, 1)) |
|
461 return Matrix (); |
|
462 |
1
|
463 int ar = a.rows (); |
|
464 int ac = a.columns (); |
|
465 |
|
466 Matrix c (ar, ac); |
|
467 |
|
468 for (int j = 0; j < ac; j++) |
|
469 for (int i = 0; i < ar; i++) |
|
470 { |
|
471 switch (op) |
|
472 { |
|
473 case Matrix_LT: |
|
474 c.elem (i, j) = a.elem (i, j) < real (b.elem (i, j)); |
|
475 break; |
|
476 case Matrix_LE: |
|
477 c.elem (i, j) = a.elem (i, j) <= real (b.elem (i, j)); |
|
478 break; |
|
479 case Matrix_EQ: |
|
480 c.elem (i, j) = a.elem (i, j) == b.elem (i, j); |
|
481 break; |
|
482 case Matrix_GE: |
|
483 c.elem (i, j) = a.elem (i, j) >= real (b.elem (i, j)); |
|
484 break; |
|
485 case Matrix_GT: |
|
486 c.elem (i, j) = a.elem (i, j) > real (b.elem (i, j)); |
|
487 break; |
|
488 case Matrix_NE: |
|
489 c.elem (i, j) = a.elem (i, j) != b.elem (i, j); |
|
490 break; |
|
491 case Matrix_AND: |
|
492 c.elem (i, j) = a.elem (i, j) && (b.elem (i, j) != 0.0); |
|
493 break; |
|
494 case Matrix_OR: |
|
495 c.elem (i, j) = a.elem (i, j) || (b.elem (i, j) != 0.0); |
|
496 break; |
|
497 default: |
|
498 panic_impossible (); |
|
499 break; |
|
500 } |
|
501 } |
|
502 return c; |
|
503 } |
|
504 |
|
505 /* 7 */ |
|
506 static Matrix |
164
|
507 mx_stupid_bool_op (Matrix_bool_op op, const Complex& s, const Matrix& a) |
1
|
508 { |
|
509 int ar = a.rows (); |
|
510 int ac = a.columns (); |
143
|
511 |
1
|
512 Matrix t (ar, ac); |
143
|
513 |
1
|
514 for (int j = 0; j < ac; j++) |
|
515 for (int i = 0; i < ar; i++) |
|
516 { |
|
517 switch (op) |
|
518 { |
|
519 case Matrix_LT: |
|
520 t.elem (i,j) = real (s) < a.elem (i,j); |
|
521 break; |
|
522 case Matrix_LE: |
|
523 t.elem (i,j) = real (s) <= a.elem (i,j); |
|
524 break; |
|
525 case Matrix_EQ: |
|
526 t.elem (i,j) = s == a.elem (i,j); |
|
527 break; |
|
528 case Matrix_GE: |
|
529 t.elem (i,j) = real (s) >= a.elem (i,j); |
|
530 break; |
|
531 case Matrix_GT: |
|
532 t.elem (i,j) = real (s) > a.elem (i,j); |
|
533 break; |
|
534 case Matrix_NE: |
|
535 t.elem (i,j) = s != a.elem (i,j); |
|
536 break; |
|
537 case Matrix_AND: |
|
538 t.elem (i,j) = (s != 0.0) && a.elem (i,j); |
|
539 break; |
|
540 case Matrix_OR: |
|
541 t.elem (i,j) = (s != 0.0) || a.elem (i,j); |
|
542 break; |
|
543 default: |
|
544 panic_impossible (); |
|
545 break; |
|
546 } |
|
547 } |
143
|
548 |
1
|
549 return t; |
|
550 } |
|
551 |
|
552 /* 8 */ |
|
553 static Matrix |
164
|
554 mx_stupid_bool_op (Matrix_bool_op op, const Complex& s, const ComplexMatrix& a) |
1
|
555 { |
|
556 int ar = a.rows (); |
|
557 int ac = a.columns (); |
143
|
558 |
1
|
559 Matrix t (ar, ac); |
143
|
560 |
1
|
561 for (int j = 0; j < ac; j++) |
|
562 for (int i = 0; i < ar; i++) |
|
563 { |
|
564 switch (op) |
|
565 { |
|
566 case Matrix_LT: |
|
567 t.elem (i,j) = real (s) < real (a.elem (i,j)); |
|
568 break; |
|
569 case Matrix_LE: |
|
570 t.elem (i,j) = real (s) <= real (a.elem (i,j)); |
|
571 break; |
|
572 case Matrix_EQ: |
|
573 t.elem (i,j) = s == a.elem (i,j); |
|
574 break; |
|
575 case Matrix_GE: |
|
576 t.elem (i,j) = real (s) >= real (a.elem (i,j)); |
|
577 break; |
|
578 case Matrix_GT: |
|
579 t.elem (i,j) = real (s) > real (a.elem (i,j)); |
|
580 break; |
|
581 case Matrix_NE: |
|
582 t.elem (i,j) = s != a.elem (i,j); |
|
583 break; |
|
584 case Matrix_AND: |
|
585 t.elem (i,j) = (s != 0.0) && (a.elem (i,j) != 0.0); |
|
586 break; |
|
587 case Matrix_OR: |
|
588 t.elem (i,j) = (s != 0.0) || (a.elem (i,j) != 0.0); |
|
589 break; |
|
590 default: |
|
591 panic_impossible (); |
|
592 break; |
|
593 } |
|
594 } |
143
|
595 |
1
|
596 return t; |
|
597 } |
|
598 |
|
599 /* 9 */ |
|
600 static Matrix |
164
|
601 mx_stupid_bool_op (Matrix_bool_op op, const ComplexMatrix& a, double s) |
1
|
602 { |
|
603 int ar = a.rows (); |
|
604 int ac = a.columns (); |
143
|
605 |
1
|
606 Matrix t (ar, ac); |
143
|
607 |
1
|
608 for (int j = 0; j < ac; j++) |
|
609 for (int i = 0; i < ar; i++) |
|
610 { |
|
611 switch (op) |
|
612 { |
|
613 case Matrix_LT: |
|
614 t.elem (i,j) = real (a.elem (i,j)) < s; |
|
615 break; |
|
616 case Matrix_LE: |
|
617 t.elem (i,j) = real (a.elem (i,j)) <= s; |
|
618 break; |
|
619 case Matrix_EQ: |
|
620 t.elem (i,j) = a.elem (i,j) == s; |
|
621 break; |
|
622 case Matrix_GE: |
|
623 t.elem (i,j) = real (a.elem (i,j)) >= s; |
|
624 break; |
|
625 case Matrix_GT: |
|
626 t.elem (i,j) = real (a.elem (i,j)) > s; |
|
627 break; |
|
628 case Matrix_NE: |
|
629 t.elem (i,j) = a.elem (i,j) != s; |
|
630 break; |
|
631 case Matrix_AND: |
|
632 t.elem (i,j) = (a.elem (i,j) != 0.0) && s; |
|
633 break; |
|
634 case Matrix_OR: |
|
635 t.elem (i,j) = (a.elem (i,j) != 0.0) || s; |
|
636 break; |
|
637 default: |
|
638 panic_impossible (); |
|
639 break; |
|
640 } |
|
641 } |
143
|
642 |
1
|
643 return t; |
|
644 } |
|
645 |
|
646 /* 10 */ |
|
647 static Matrix |
164
|
648 mx_stupid_bool_op (Matrix_bool_op op, const ComplexMatrix& a, const Complex& s) |
1
|
649 { |
|
650 int ar = a.rows (); |
|
651 int ac = a.columns (); |
143
|
652 |
1
|
653 Matrix t (ar, ac); |
143
|
654 |
1
|
655 for (int j = 0; j < ac; j++) |
|
656 for (int i = 0; i < ar; i++) |
|
657 { |
|
658 switch (op) |
|
659 { |
|
660 case Matrix_LT: |
|
661 t.elem (i,j) = real (a.elem (i,j)) < real (s); |
|
662 break; |
|
663 case Matrix_LE: |
|
664 t.elem (i,j) = real (a.elem (i,j)) <= real (s); |
|
665 break; |
|
666 case Matrix_EQ: |
|
667 t.elem (i,j) = a.elem (i,j) == s; |
|
668 break; |
|
669 case Matrix_GE: |
|
670 t.elem (i,j) = real (a.elem (i,j)) >= real (s); |
|
671 break; |
|
672 case Matrix_GT: |
|
673 t.elem (i,j) = real (a.elem (i,j)) > real (s); |
|
674 break; |
|
675 case Matrix_NE: |
|
676 t.elem (i,j) = a.elem (i,j) != s; |
|
677 break; |
|
678 case Matrix_AND: |
|
679 t.elem (i,j) = (a.elem (i,j) != 0.0) && (s != 0.0); |
|
680 break; |
|
681 case Matrix_OR: |
|
682 t.elem (i,j) = (a.elem (i,j) != 0.0) || (s != 0.0); |
|
683 break; |
|
684 default: |
|
685 panic_impossible (); |
|
686 break; |
|
687 } |
|
688 } |
143
|
689 |
1
|
690 return t; |
|
691 } |
|
692 |
|
693 /* 11 */ |
|
694 static Matrix |
164
|
695 mx_stupid_bool_op (Matrix_bool_op op, const ComplexMatrix& a, const Matrix& b) |
1
|
696 { |
143
|
697 if (! m_add_conform (a, b, 1)) |
|
698 return Matrix (); |
|
699 |
1
|
700 int ar = a.rows (); |
|
701 int ac = a.columns (); |
|
702 |
|
703 Matrix c (ar, ac); |
|
704 |
|
705 for (int j = 0; j < ac; j++) |
|
706 for (int i = 0; i < ar; i++) |
|
707 { |
|
708 switch (op) |
|
709 { |
|
710 case Matrix_LT: |
|
711 c.elem (i, j) = real (a.elem (i, j)) < b.elem (i, j); |
|
712 break; |
|
713 case Matrix_LE: |
|
714 c.elem (i, j) = real (a.elem (i, j)) <= b.elem (i, j); |
|
715 break; |
|
716 case Matrix_EQ: |
|
717 c.elem (i, j) = a.elem (i, j) == b.elem (i, j); |
|
718 break; |
|
719 case Matrix_GE: |
|
720 c.elem (i, j) = real (a.elem (i, j)) >= b.elem (i, j); |
|
721 break; |
|
722 case Matrix_GT: |
|
723 c.elem (i, j) = real (a.elem (i, j)) > b.elem (i, j); |
|
724 break; |
|
725 case Matrix_NE: |
|
726 c.elem (i, j) = a.elem (i, j) != b.elem (i, j); |
|
727 break; |
|
728 case Matrix_AND: |
|
729 c.elem (i, j) = (a.elem (i, j) != 0.0) && b.elem (i, j); |
|
730 break; |
|
731 case Matrix_OR: |
|
732 c.elem (i, j) = (a.elem (i, j) != 0.0) || b.elem (i, j); |
|
733 break; |
|
734 default: |
|
735 panic_impossible (); |
|
736 break; |
|
737 } |
|
738 } |
|
739 return c; |
|
740 } |
|
741 |
|
742 /* 12 */ |
|
743 static Matrix |
164
|
744 mx_stupid_bool_op (Matrix_bool_op op, const ComplexMatrix& a, |
|
745 const ComplexMatrix& b) |
1
|
746 { |
143
|
747 if (! m_add_conform (a, b, 1)) |
|
748 return Matrix (); |
|
749 |
1
|
750 int ar = a.rows (); |
|
751 int ac = a.columns (); |
|
752 |
|
753 Matrix c (ar, ac); |
|
754 |
|
755 for (int j = 0; j < ac; j++) |
|
756 for (int i = 0; i < ar; i++) |
|
757 { |
|
758 switch (op) |
|
759 { |
|
760 case Matrix_LT: |
|
761 c.elem (i, j) = real (a.elem (i, j)) < real (b.elem (i, j)); |
|
762 break; |
|
763 case Matrix_LE: |
|
764 c.elem (i, j) = real (a.elem (i, j)) <= real (b.elem (i, j)); |
|
765 break; |
|
766 case Matrix_EQ: |
|
767 c.elem (i, j) = a.elem (i, j) == b.elem (i, j); |
|
768 break; |
|
769 case Matrix_GE: |
|
770 c.elem (i, j) = real (a.elem (i, j)) >= real (b.elem (i, j)); |
|
771 break; |
|
772 case Matrix_GT: |
|
773 c.elem (i, j) = real (a.elem (i, j)) > real (b.elem (i, j)); |
|
774 break; |
|
775 case Matrix_NE: |
|
776 c.elem (i, j) = a.elem (i, j) != b.elem (i, j); |
|
777 break; |
|
778 case Matrix_AND: |
|
779 c.elem (i, j) = (a.elem (i, j) != 0.0) && (b.elem (i, j) != 0.0); |
|
780 break; |
|
781 case Matrix_OR: |
|
782 c.elem (i, j) = (a.elem (i, j) != 0.0) || (b.elem (i, j) != 0.0); |
|
783 break; |
|
784 default: |
|
785 panic_impossible (); |
|
786 break; |
|
787 } |
|
788 } |
143
|
789 |
1
|
790 return c; |
|
791 } |
|
792 |
|
793 /* |
|
794 * Unary operations. One for each numeric data type: |
|
795 * |
|
796 * scalar |
|
797 * complex_scalar |
|
798 * matrix |
|
799 * complex_matrix |
|
800 * |
|
801 */ |
|
802 |
|
803 tree_constant |
|
804 do_unary_op (double d, tree::expression_type t) |
|
805 { |
|
806 double result = 0.0; |
143
|
807 |
1
|
808 switch (t) |
|
809 { |
|
810 case tree::not: |
|
811 result = (! d); |
|
812 break; |
|
813 case tree::uminus: |
|
814 result = -d; |
|
815 break; |
|
816 case tree::hermitian: |
|
817 case tree::transpose: |
|
818 result = d; |
|
819 break; |
|
820 default: |
|
821 panic_impossible (); |
|
822 break; |
|
823 } |
|
824 |
|
825 return tree_constant (result); |
|
826 } |
|
827 |
|
828 tree_constant |
164
|
829 do_unary_op (const Matrix& a, tree::expression_type t) |
1
|
830 { |
|
831 Matrix result; |
143
|
832 |
1
|
833 switch (t) |
|
834 { |
|
835 case tree::not: |
|
836 result = (! a); |
|
837 break; |
|
838 case tree::uminus: |
|
839 result = -a; |
|
840 break; |
|
841 case tree::hermitian: |
|
842 case tree::transpose: |
|
843 result = a.transpose (); |
|
844 break; |
|
845 default: |
|
846 panic_impossible (); |
|
847 break; |
|
848 } |
|
849 |
|
850 return tree_constant (result); |
|
851 } |
|
852 |
|
853 tree_constant |
164
|
854 do_unary_op (const Complex& c, tree::expression_type t) |
1
|
855 { |
|
856 Complex result = 0.0; |
143
|
857 |
1
|
858 switch (t) |
|
859 { |
|
860 case tree::not: |
|
861 result = (c == 0.0); |
|
862 break; |
|
863 case tree::uminus: |
|
864 result = -c; |
|
865 break; |
|
866 case tree::hermitian: |
|
867 result = conj (c); |
|
868 break; |
|
869 case tree::transpose: |
|
870 result = c; |
|
871 break; |
|
872 default: |
|
873 panic_impossible (); |
|
874 break; |
|
875 } |
|
876 |
|
877 return tree_constant (result); |
|
878 } |
|
879 |
|
880 tree_constant |
164
|
881 do_unary_op (const ComplexMatrix& a, tree::expression_type t) |
1
|
882 { |
|
883 ComplexMatrix result; |
143
|
884 |
1
|
885 switch (t) |
|
886 { |
|
887 case tree::not: |
|
888 result = (! a); |
|
889 break; |
|
890 case tree::uminus: |
|
891 result = -a; |
|
892 break; |
|
893 case tree::hermitian: |
|
894 result = a.hermitian (); |
|
895 break; |
|
896 case tree::transpose: |
|
897 result = a.transpose (); |
|
898 break; |
|
899 default: |
|
900 panic_impossible (); |
|
901 break; |
|
902 } |
|
903 |
|
904 return tree_constant (result); |
|
905 } |
|
906 |
|
907 /* |
|
908 * Binary operations. One for each type combination, in the order |
|
909 * given here: |
|
910 * |
|
911 * op2 \ op1: s m cs cm |
|
912 * +-- +---+---+----+----+ |
|
913 * scalar | | 1 | 5 | 9 | 13 | |
|
914 * +---+---+----+----+ |
|
915 * matrix | 2 | 6 | 10 | 14 | |
|
916 * +---+---+----+----+ |
|
917 * complex_scalar | 3 | 7 | 11 | 15 | |
|
918 * +---+---+----+----+ |
|
919 * complex_matrix | 4 | 8 | 12 | 16 | |
|
920 * +---+---+----+----+ |
|
921 */ |
|
922 |
|
923 /* 1 */ |
|
924 tree_constant |
|
925 do_binary_op (double a, double b, tree::expression_type t) |
|
926 { |
|
927 double result = 0.0; |
143
|
928 |
1
|
929 switch (t) |
|
930 { |
|
931 case tree::add: |
|
932 result = a + b; |
|
933 break; |
|
934 case tree::subtract: |
|
935 result = a - b; |
|
936 break; |
|
937 case tree::multiply: |
|
938 case tree::el_mul: |
|
939 result = a * b; |
|
940 break; |
|
941 case tree::divide: |
|
942 case tree::el_div: |
|
943 if (b == 0.0) |
|
944 DIVIDE_BY_ZERO_ERROR; |
|
945 result = a / b; |
|
946 break; |
|
947 case tree::leftdiv: |
|
948 case tree::el_leftdiv: |
|
949 if (a == 0.0) |
|
950 DIVIDE_BY_ZERO_ERROR; |
|
951 result = b / a; |
|
952 break; |
|
953 case tree::power: |
|
954 case tree::elem_pow: |
|
955 return xpow (a, b); |
|
956 break; |
|
957 case tree::cmp_lt: |
|
958 result = a < b; |
|
959 break; |
|
960 case tree::cmp_le: |
|
961 result = a <= b; |
|
962 break; |
|
963 case tree::cmp_eq: |
|
964 result = a == b; |
|
965 break; |
|
966 case tree::cmp_ge: |
|
967 result = a >= b; |
|
968 break; |
|
969 case tree::cmp_gt: |
|
970 result = a > b; |
|
971 break; |
|
972 case tree::cmp_ne: |
|
973 result = a != b; |
|
974 break; |
|
975 case tree::and: |
|
976 result = (a && b); |
|
977 break; |
|
978 case tree::or: |
|
979 result = (a || b); |
|
980 break; |
|
981 default: |
|
982 panic_impossible (); |
|
983 break; |
|
984 } |
|
985 |
143
|
986 if (error_state) |
|
987 return tree_constant (); |
|
988 |
1
|
989 return tree_constant (result); |
|
990 } |
|
991 |
|
992 /* 2 */ |
|
993 tree_constant |
164
|
994 do_binary_op (double a, const Matrix& b, tree::expression_type t) |
1
|
995 { |
|
996 Matrix result; |
143
|
997 |
1
|
998 switch (t) |
|
999 { |
|
1000 case tree::add: |
|
1001 result = a + b; |
|
1002 break; |
|
1003 case tree::subtract: |
|
1004 result = a - b; |
|
1005 break; |
|
1006 case tree::el_leftdiv: |
|
1007 case tree::leftdiv: |
|
1008 if (a == 0.0) |
|
1009 DIVIDE_BY_ZERO_ERROR; |
|
1010 a = 1.0 / a; |
|
1011 // fall through... |
|
1012 case tree::multiply: |
|
1013 case tree::el_mul: |
|
1014 result = a * b; |
|
1015 break; |
|
1016 case tree::el_div: |
|
1017 return x_el_div (a, b); |
|
1018 break; |
|
1019 case tree::divide: |
143
|
1020 gripe_nonconformant (1, 1, b.rows (), b.columns ()); |
1
|
1021 break; |
|
1022 case tree::power: |
|
1023 return xpow (a, b); |
|
1024 break; |
|
1025 case tree::elem_pow: |
|
1026 return elem_xpow (a, b); |
|
1027 break; |
|
1028 case tree::cmp_lt: |
|
1029 result = mx_stupid_bool_op (Matrix_LT, a, b); |
|
1030 break; |
|
1031 case tree::cmp_le: |
|
1032 result = mx_stupid_bool_op (Matrix_LE, a, b); |
|
1033 break; |
|
1034 case tree::cmp_eq: |
|
1035 result = mx_stupid_bool_op (Matrix_EQ, a, b); |
|
1036 break; |
|
1037 case tree::cmp_ge: |
|
1038 result = mx_stupid_bool_op (Matrix_GE, a, b); |
|
1039 break; |
|
1040 case tree::cmp_gt: |
|
1041 result = mx_stupid_bool_op (Matrix_GT, a, b); |
|
1042 break; |
|
1043 case tree::cmp_ne: |
|
1044 result = mx_stupid_bool_op (Matrix_NE, a, b); |
|
1045 break; |
|
1046 case tree::and: |
|
1047 result = mx_stupid_bool_op (Matrix_AND, a, b); |
|
1048 break; |
|
1049 case tree::or: |
|
1050 result = mx_stupid_bool_op (Matrix_OR, a, b); |
|
1051 break; |
|
1052 default: |
|
1053 panic_impossible (); |
|
1054 break; |
|
1055 } |
|
1056 |
143
|
1057 if (error_state) |
|
1058 return tree_constant (); |
|
1059 |
1
|
1060 return tree_constant (result); |
|
1061 } |
|
1062 |
|
1063 /* 3 */ |
|
1064 tree_constant |
164
|
1065 do_binary_op (double a, const Complex& b, tree::expression_type t) |
1
|
1066 { |
|
1067 enum RT { RT_unknown, RT_real, RT_complex }; |
|
1068 RT result_type = RT_unknown; |
|
1069 |
|
1070 double result = 0.0; |
|
1071 Complex complex_result; |
143
|
1072 |
1
|
1073 switch (t) |
|
1074 { |
|
1075 case tree::add: |
|
1076 result_type = RT_complex; |
|
1077 complex_result = a + b; |
|
1078 break; |
|
1079 case tree::subtract: |
|
1080 result_type = RT_complex; |
|
1081 complex_result = a - b; |
|
1082 break; |
|
1083 case tree::multiply: |
|
1084 case tree::el_mul: |
|
1085 result_type = RT_complex; |
|
1086 complex_result = a * b; |
|
1087 break; |
|
1088 case tree::divide: |
|
1089 case tree::el_div: |
|
1090 result_type = RT_complex; |
|
1091 if (b == 0.0) |
|
1092 DIVIDE_BY_ZERO_ERROR; |
|
1093 complex_result = a / b; |
|
1094 break; |
|
1095 case tree::leftdiv: |
|
1096 case tree::el_leftdiv: |
|
1097 result_type = RT_complex; |
|
1098 if (a == 0.0) |
|
1099 DIVIDE_BY_ZERO_ERROR; |
|
1100 complex_result = b / a; |
|
1101 break; |
|
1102 case tree::power: |
|
1103 case tree::elem_pow: |
|
1104 return xpow (a, b); |
|
1105 break; |
|
1106 case tree::cmp_lt: |
|
1107 result_type = RT_real; |
|
1108 result = a < real (b); |
|
1109 break; |
|
1110 case tree::cmp_le: |
|
1111 result_type = RT_real; |
|
1112 result = a <= real (b); |
|
1113 break; |
|
1114 case tree::cmp_eq: |
|
1115 result_type = RT_real; |
|
1116 result = a == b; |
|
1117 break; |
|
1118 case tree::cmp_ge: |
|
1119 result_type = RT_real; |
|
1120 result = a >= real (b); |
|
1121 break; |
|
1122 case tree::cmp_gt: |
|
1123 result_type = RT_real; |
|
1124 result = a > real (b); |
|
1125 break; |
|
1126 case tree::cmp_ne: |
|
1127 result_type = RT_real; |
|
1128 result = a != b; |
|
1129 break; |
|
1130 case tree::and: |
|
1131 result_type = RT_real; |
|
1132 result = (a && (b != 0.0)); |
|
1133 break; |
|
1134 case tree::or: |
|
1135 result_type = RT_real; |
|
1136 result = (a || (b != 0.0)); |
|
1137 break; |
|
1138 default: |
|
1139 panic_impossible (); |
|
1140 break; |
|
1141 } |
|
1142 |
143
|
1143 if (error_state) |
|
1144 return tree_constant (); |
|
1145 |
1
|
1146 assert (result_type != RT_unknown); |
143
|
1147 |
1
|
1148 if (result_type == RT_real) |
|
1149 return tree_constant (result); |
|
1150 else |
|
1151 return tree_constant (complex_result); |
|
1152 } |
|
1153 |
|
1154 /* 4 */ |
|
1155 tree_constant |
164
|
1156 do_binary_op (double a, const ComplexMatrix& b, tree::expression_type t) |
1
|
1157 { |
|
1158 enum RT { RT_unknown, RT_real, RT_complex }; |
|
1159 RT result_type = RT_unknown; |
|
1160 |
|
1161 Matrix result; |
|
1162 ComplexMatrix complex_result; |
143
|
1163 |
1
|
1164 switch (t) |
|
1165 { |
|
1166 case tree::add: |
|
1167 result_type = RT_complex; |
|
1168 complex_result = a + b; |
|
1169 break; |
|
1170 case tree::subtract: |
|
1171 result_type = RT_complex; |
|
1172 complex_result = a - b; |
|
1173 break; |
|
1174 case tree::el_leftdiv: |
|
1175 case tree::leftdiv: |
|
1176 if (a == 0.0) |
|
1177 DIVIDE_BY_ZERO_ERROR; |
|
1178 a = 1.0 / a; |
|
1179 // fall through... |
|
1180 case tree::multiply: |
|
1181 case tree::el_mul: |
|
1182 result_type = RT_complex; |
|
1183 complex_result = a * b; |
|
1184 break; |
|
1185 case tree::el_div: |
|
1186 return x_el_div (a, b); |
|
1187 break; |
|
1188 case tree::divide: |
143
|
1189 gripe_nonconformant (1, 1, b.rows (), b.columns ()); |
1
|
1190 break; |
|
1191 case tree::power: |
|
1192 return xpow (a, b); |
|
1193 break; |
|
1194 case tree::elem_pow: |
|
1195 return elem_xpow (a, b); |
|
1196 break; |
|
1197 case tree::cmp_lt: |
|
1198 result_type = RT_real; |
|
1199 result = mx_stupid_bool_op (Matrix_LT, a, b); |
|
1200 break; |
|
1201 case tree::cmp_le: |
|
1202 result_type = RT_real; |
|
1203 result = mx_stupid_bool_op (Matrix_LE, a, b); |
|
1204 break; |
|
1205 case tree::cmp_eq: |
|
1206 result_type = RT_real; |
|
1207 result = mx_stupid_bool_op (Matrix_EQ, a, b); |
|
1208 break; |
|
1209 case tree::cmp_ge: |
|
1210 result_type = RT_real; |
|
1211 result = mx_stupid_bool_op (Matrix_GE, a, b); |
|
1212 break; |
|
1213 case tree::cmp_gt: |
|
1214 result_type = RT_real; |
|
1215 result = mx_stupid_bool_op (Matrix_GT, a, b); |
|
1216 break; |
|
1217 case tree::cmp_ne: |
|
1218 result_type = RT_real; |
|
1219 result = mx_stupid_bool_op (Matrix_NE, a, b); |
|
1220 break; |
|
1221 case tree::and: |
|
1222 result_type = RT_real; |
|
1223 result = mx_stupid_bool_op (Matrix_AND, a, b); |
|
1224 break; |
|
1225 case tree::or: |
|
1226 result_type = RT_real; |
|
1227 result = mx_stupid_bool_op (Matrix_OR, a, b); |
|
1228 break; |
|
1229 default: |
|
1230 panic_impossible (); |
|
1231 break; |
|
1232 } |
|
1233 |
143
|
1234 if (error_state) |
|
1235 return tree_constant (); |
|
1236 |
1
|
1237 assert (result_type != RT_unknown); |
143
|
1238 |
1
|
1239 if (result_type == RT_real) |
|
1240 return tree_constant (result); |
|
1241 else |
|
1242 return tree_constant (complex_result); |
|
1243 } |
|
1244 |
|
1245 /* 5 */ |
|
1246 tree_constant |
164
|
1247 do_binary_op (const Matrix& a, double b, tree::expression_type t) |
1
|
1248 { |
|
1249 Matrix result; |
143
|
1250 |
1
|
1251 switch (t) |
|
1252 { |
|
1253 case tree::add: |
|
1254 result = a + b; |
|
1255 break; |
|
1256 case tree::subtract: |
|
1257 result = a - b; |
|
1258 break; |
|
1259 case tree::multiply: |
|
1260 case tree::el_mul: |
|
1261 result = a * b; |
|
1262 break; |
|
1263 case tree::divide: |
|
1264 case tree::el_div: |
|
1265 result = a / b; |
|
1266 break; |
|
1267 case tree::el_leftdiv: |
|
1268 return x_el_div (b, a); |
|
1269 break; |
|
1270 case tree::leftdiv: |
143
|
1271 gripe_nonconformant (a.rows (), a.columns (), 1, 1); |
1
|
1272 break; |
|
1273 case tree::power: |
|
1274 return xpow (a, b); |
|
1275 break; |
|
1276 case tree::elem_pow: |
|
1277 return elem_xpow (a, b); |
|
1278 break; |
|
1279 case tree::cmp_lt: |
|
1280 result = mx_stupid_bool_op (Matrix_LT, a, b); |
|
1281 break; |
|
1282 case tree::cmp_le: |
|
1283 result = mx_stupid_bool_op (Matrix_LE, a, b); |
|
1284 break; |
|
1285 case tree::cmp_eq: |
|
1286 result = mx_stupid_bool_op (Matrix_EQ, a, b); |
|
1287 break; |
|
1288 case tree::cmp_ge: |
|
1289 result = mx_stupid_bool_op (Matrix_GE, a, b); |
|
1290 break; |
|
1291 case tree::cmp_gt: |
|
1292 result = mx_stupid_bool_op (Matrix_GT, a, b); |
|
1293 break; |
|
1294 case tree::cmp_ne: |
|
1295 result = mx_stupid_bool_op (Matrix_NE, a, b); |
|
1296 break; |
|
1297 case tree::and: |
|
1298 result = mx_stupid_bool_op (Matrix_AND, a, b); |
|
1299 break; |
|
1300 case tree::or: |
|
1301 result = mx_stupid_bool_op (Matrix_OR, a, b); |
|
1302 break; |
|
1303 default: |
|
1304 panic_impossible (); |
|
1305 break; |
|
1306 } |
|
1307 |
143
|
1308 if (error_state) |
|
1309 return tree_constant (); |
|
1310 |
1
|
1311 return tree_constant (result); |
|
1312 } |
|
1313 |
|
1314 /* 6 */ |
|
1315 tree_constant |
164
|
1316 do_binary_op (const Matrix& a, const Matrix& b, tree::expression_type t) |
1
|
1317 { |
|
1318 Matrix result; |
|
1319 |
|
1320 switch (t) |
|
1321 { |
|
1322 case tree::add: |
|
1323 if (m_add_conform (a, b, 1)) |
|
1324 result = a + b; |
|
1325 break; |
|
1326 case tree::subtract: |
|
1327 if (m_add_conform (a, b, 1)) |
|
1328 result = a - b; |
|
1329 break; |
|
1330 case tree::el_mul: |
|
1331 if (m_add_conform (a, b, 1)) |
240
|
1332 result = product (a, b); |
1
|
1333 break; |
|
1334 case tree::multiply: |
|
1335 if (m_mul_conform (a, b, 1)) |
|
1336 result = a * b; |
|
1337 break; |
|
1338 case tree::el_div: |
|
1339 if (m_add_conform (a, b, 1)) |
240
|
1340 result = quotient (a, b); |
1
|
1341 break; |
|
1342 case tree::el_leftdiv: |
|
1343 if (m_add_conform (a, b, 1)) |
240
|
1344 result = quotient (b, a); |
1
|
1345 break; |
|
1346 case tree::leftdiv: |
|
1347 return xleftdiv (a, b); |
|
1348 break; |
|
1349 case tree::divide: |
|
1350 return xdiv (a, b); |
|
1351 break; |
|
1352 case tree::power: |
|
1353 error ("can't do A ^ B for A and B both matrices"); |
|
1354 break; |
|
1355 case tree::elem_pow: |
|
1356 if (m_add_conform (a, b, 1)) |
|
1357 return elem_xpow (a, b); |
|
1358 break; |
|
1359 case tree::cmp_lt: |
|
1360 if (m_add_conform (a, b, 1)) |
|
1361 result = mx_stupid_bool_op (Matrix_LT, a, b); |
|
1362 break; |
|
1363 case tree::cmp_le: |
|
1364 if (m_add_conform (a, b, 1)) |
|
1365 result = mx_stupid_bool_op (Matrix_LE, a, b); |
|
1366 break; |
|
1367 case tree::cmp_eq: |
|
1368 if (m_add_conform (a, b, 1)) |
|
1369 result = mx_stupid_bool_op (Matrix_EQ, a, b); |
|
1370 break; |
|
1371 case tree::cmp_ge: |
|
1372 if (m_add_conform (a, b, 1)) |
|
1373 result = mx_stupid_bool_op (Matrix_GE, a, b); |
|
1374 break; |
|
1375 case tree::cmp_gt: |
|
1376 if (m_add_conform (a, b, 1)) |
|
1377 result = mx_stupid_bool_op (Matrix_GT, a, b); |
|
1378 break; |
|
1379 case tree::cmp_ne: |
|
1380 if (m_add_conform (a, b, 1)) |
|
1381 result = mx_stupid_bool_op (Matrix_NE, a, b); |
|
1382 break; |
|
1383 case tree::and: |
|
1384 if (m_add_conform (a, b, 1)) |
|
1385 result = mx_stupid_bool_op (Matrix_AND, a, b); |
|
1386 break; |
|
1387 case tree::or: |
|
1388 if (m_add_conform (a, b, 1)) |
|
1389 result = mx_stupid_bool_op (Matrix_OR, a, b); |
|
1390 break; |
|
1391 default: |
|
1392 panic_impossible (); |
|
1393 break; |
|
1394 } |
|
1395 |
143
|
1396 if (error_state) |
1
|
1397 return tree_constant (); |
143
|
1398 |
|
1399 return tree_constant (result); |
1
|
1400 } |
|
1401 |
|
1402 /* 7 */ |
|
1403 tree_constant |
164
|
1404 do_binary_op (const Matrix& a, const Complex& b, tree::expression_type t) |
1
|
1405 { |
|
1406 enum RT { RT_unknown, RT_real, RT_complex }; |
|
1407 RT result_type = RT_unknown; |
|
1408 |
|
1409 Matrix result; |
|
1410 ComplexMatrix complex_result; |
143
|
1411 |
1
|
1412 switch (t) |
|
1413 { |
|
1414 case tree::add: |
|
1415 result_type = RT_complex; |
|
1416 complex_result = a + b; |
|
1417 break; |
|
1418 case tree::subtract: |
|
1419 result_type = RT_complex; |
|
1420 complex_result = a - b; |
|
1421 break; |
|
1422 case tree::multiply: |
|
1423 case tree::el_mul: |
|
1424 result_type = RT_complex; |
|
1425 complex_result = a * b; |
|
1426 break; |
|
1427 case tree::divide: |
|
1428 case tree::el_div: |
|
1429 result_type = RT_complex; |
|
1430 complex_result = a / b; |
|
1431 break; |
|
1432 case tree::el_leftdiv: |
|
1433 return x_el_div (b, a); |
|
1434 break; |
|
1435 case tree::leftdiv: |
143
|
1436 gripe_nonconformant (a.rows (), a.columns (), 1, 1); |
1
|
1437 break; |
|
1438 case tree::power: |
|
1439 return xpow (a, b); |
|
1440 break; |
|
1441 case tree::elem_pow: |
|
1442 return elem_xpow (a, b); |
|
1443 break; |
|
1444 case tree::cmp_lt: |
|
1445 result_type = RT_real; |
|
1446 result = mx_stupid_bool_op (Matrix_LT, a, b); |
|
1447 break; |
|
1448 case tree::cmp_le: |
|
1449 result_type = RT_real; |
|
1450 result = mx_stupid_bool_op (Matrix_LE, a, b); |
|
1451 break; |
|
1452 case tree::cmp_eq: |
|
1453 result_type = RT_real; |
|
1454 result = mx_stupid_bool_op (Matrix_EQ, a, b); |
|
1455 break; |
|
1456 case tree::cmp_ge: |
|
1457 result_type = RT_real; |
|
1458 result = mx_stupid_bool_op (Matrix_GE, a, b); |
|
1459 break; |
|
1460 case tree::cmp_gt: |
|
1461 result_type = RT_real; |
|
1462 result = mx_stupid_bool_op (Matrix_GT, a, b); |
|
1463 break; |
|
1464 case tree::cmp_ne: |
|
1465 result_type = RT_real; |
|
1466 result = mx_stupid_bool_op (Matrix_NE, a, b); |
|
1467 break; |
|
1468 case tree::and: |
|
1469 result_type = RT_real; |
|
1470 result = mx_stupid_bool_op (Matrix_AND, a, b); |
|
1471 break; |
|
1472 case tree::or: |
|
1473 result_type = RT_real; |
|
1474 result = mx_stupid_bool_op (Matrix_OR, a, b); |
|
1475 break; |
|
1476 default: |
|
1477 panic_impossible (); |
|
1478 break; |
|
1479 } |
|
1480 |
143
|
1481 if (error_state) |
|
1482 return tree_constant (); |
|
1483 |
1
|
1484 assert (result_type != RT_unknown); |
143
|
1485 |
1
|
1486 if (result_type == RT_real) |
|
1487 return tree_constant (result); |
|
1488 else |
|
1489 return tree_constant (complex_result); |
|
1490 } |
|
1491 |
|
1492 /* 8 */ |
|
1493 tree_constant |
164
|
1494 do_binary_op (const Matrix& a, const ComplexMatrix& b, tree::expression_type t) |
1
|
1495 { |
|
1496 enum RT { RT_unknown, RT_real, RT_complex }; |
|
1497 RT result_type = RT_unknown; |
|
1498 |
|
1499 Matrix result; |
|
1500 ComplexMatrix complex_result; |
143
|
1501 |
1
|
1502 switch (t) |
|
1503 { |
|
1504 case tree::add: |
|
1505 result_type = RT_complex; |
|
1506 if (m_add_conform (a, b, 1)) |
|
1507 complex_result = a + b; |
|
1508 break; |
|
1509 case tree::subtract: |
|
1510 result_type = RT_complex; |
|
1511 if (m_add_conform (a, b, 1)) |
|
1512 complex_result = a - b; |
|
1513 break; |
|
1514 case tree::el_mul: |
|
1515 result_type = RT_complex; |
|
1516 if (m_add_conform (a, b, 1)) |
240
|
1517 complex_result = product (a, b); |
1
|
1518 break; |
|
1519 case tree::multiply: |
|
1520 result_type = RT_complex; |
|
1521 if (m_mul_conform (a, b, 1)) |
|
1522 complex_result = a * b; |
|
1523 break; |
|
1524 case tree::el_div: |
|
1525 result_type = RT_complex; |
|
1526 if (m_add_conform (a, b, 1)) |
240
|
1527 complex_result = quotient (a, b); |
1
|
1528 break; |
|
1529 case tree::el_leftdiv: |
|
1530 result_type = RT_complex; |
|
1531 if (m_add_conform (a, b, 1)) |
240
|
1532 complex_result = quotient (b, a); |
1
|
1533 break; |
|
1534 case tree::leftdiv: |
|
1535 return xleftdiv (a, b); |
|
1536 break; |
|
1537 case tree::divide: |
|
1538 return xdiv (a, b); |
|
1539 break; |
|
1540 case tree::power: |
|
1541 error ("can't do A ^ B for A and B both matrices"); |
|
1542 break; |
|
1543 case tree::elem_pow: |
|
1544 if (m_add_conform (a, b, 1)) |
|
1545 return elem_xpow (a, b); |
|
1546 break; |
|
1547 case tree::cmp_lt: |
|
1548 result_type = RT_real; |
|
1549 if (m_add_conform (a, b, 1)) |
|
1550 result = mx_stupid_bool_op (Matrix_LT, a, b); |
|
1551 break; |
|
1552 case tree::cmp_le: |
|
1553 result_type = RT_real; |
|
1554 if (m_add_conform (a, b, 1)) |
|
1555 result = mx_stupid_bool_op (Matrix_LE, a, b); |
|
1556 break; |
|
1557 case tree::cmp_eq: |
|
1558 result_type = RT_real; |
|
1559 if (m_add_conform (a, b, 1)) |
|
1560 result = mx_stupid_bool_op (Matrix_EQ, a, b); |
|
1561 break; |
|
1562 case tree::cmp_ge: |
|
1563 result_type = RT_real; |
|
1564 if (m_add_conform (a, b, 1)) |
|
1565 result = mx_stupid_bool_op (Matrix_GE, a, b); |
|
1566 break; |
|
1567 case tree::cmp_gt: |
|
1568 result_type = RT_real; |
|
1569 if (m_add_conform (a, b, 1)) |
|
1570 result = mx_stupid_bool_op (Matrix_GT, a, b); |
|
1571 break; |
|
1572 case tree::cmp_ne: |
|
1573 result_type = RT_real; |
|
1574 if (m_add_conform (a, b, 1)) |
|
1575 result = mx_stupid_bool_op (Matrix_NE, a, b); |
|
1576 break; |
|
1577 case tree::and: |
|
1578 result_type = RT_real; |
|
1579 if (m_add_conform (a, b, 1)) |
|
1580 result = mx_stupid_bool_op (Matrix_AND, a, b); |
|
1581 break; |
|
1582 case tree::or: |
|
1583 result_type = RT_real; |
|
1584 if (m_add_conform (a, b, 1)) |
|
1585 result = mx_stupid_bool_op (Matrix_OR, a, b); |
|
1586 break; |
|
1587 default: |
|
1588 panic_impossible (); |
|
1589 break; |
|
1590 } |
|
1591 |
143
|
1592 if (error_state) |
|
1593 return tree_constant (); |
|
1594 |
1
|
1595 assert (result_type != RT_unknown); |
143
|
1596 |
1
|
1597 if (result_type == RT_real) |
|
1598 return tree_constant (result); |
|
1599 else |
|
1600 return tree_constant (complex_result); |
|
1601 } |
|
1602 |
|
1603 /* 9 */ |
|
1604 tree_constant |
164
|
1605 do_binary_op (const Complex& a, double b, tree::expression_type t) |
1
|
1606 { |
|
1607 enum RT { RT_unknown, RT_real, RT_complex }; |
|
1608 RT result_type = RT_unknown; |
|
1609 |
|
1610 double result = 0.0; |
|
1611 Complex complex_result; |
143
|
1612 |
1
|
1613 switch (t) |
|
1614 { |
|
1615 case tree::add: |
|
1616 result_type = RT_complex; |
|
1617 complex_result = a + b; |
|
1618 break; |
|
1619 case tree::subtract: |
|
1620 result_type = RT_complex; |
|
1621 complex_result = a - b; |
|
1622 break; |
|
1623 case tree::multiply: |
|
1624 case tree::el_mul: |
|
1625 result_type = RT_complex; |
|
1626 complex_result = a * b; |
|
1627 break; |
|
1628 case tree::divide: |
|
1629 case tree::el_div: |
|
1630 result_type = RT_complex; |
|
1631 if (b == 0.0) |
|
1632 DIVIDE_BY_ZERO_ERROR; |
|
1633 complex_result = a / b; |
|
1634 break; |
|
1635 case tree::leftdiv: |
|
1636 case tree::el_leftdiv: |
|
1637 result_type = RT_complex; |
|
1638 if (a == 0.0) |
|
1639 DIVIDE_BY_ZERO_ERROR; |
|
1640 complex_result = b / a; |
|
1641 break; |
|
1642 case tree::power: |
|
1643 case tree::elem_pow: |
|
1644 return xpow (a, b); |
|
1645 break; |
|
1646 case tree::cmp_lt: |
|
1647 result_type = RT_real; |
|
1648 result = real (a) < b; |
|
1649 break; |
|
1650 case tree::cmp_le: |
|
1651 result_type = RT_real; |
|
1652 result = real (a) <= b; |
|
1653 break; |
|
1654 case tree::cmp_eq: |
|
1655 result_type = RT_real; |
|
1656 result = a == b; |
|
1657 break; |
|
1658 case tree::cmp_ge: |
|
1659 result_type = RT_real; |
|
1660 result = real (a) >= b; |
|
1661 break; |
|
1662 case tree::cmp_gt: |
|
1663 result_type = RT_real; |
|
1664 result = real (a) > b; |
|
1665 break; |
|
1666 case tree::cmp_ne: |
|
1667 result_type = RT_real; |
|
1668 result = a != b; |
|
1669 break; |
|
1670 case tree::and: |
|
1671 result_type = RT_real; |
|
1672 result = ((a != 0.0) && b); |
|
1673 break; |
|
1674 case tree::or: |
|
1675 result_type = RT_real; |
|
1676 result = ((a != 0.0) || b); |
|
1677 break; |
|
1678 default: |
|
1679 panic_impossible (); |
|
1680 break; |
|
1681 } |
|
1682 |
143
|
1683 if (error_state) |
|
1684 return tree_constant (); |
|
1685 |
1
|
1686 assert (result_type != RT_unknown); |
143
|
1687 |
1
|
1688 if (result_type == RT_real) |
|
1689 return tree_constant (result); |
|
1690 else |
|
1691 return tree_constant (complex_result); |
|
1692 } |
|
1693 |
|
1694 /* 10 */ |
|
1695 tree_constant |
164
|
1696 do_binary_op (const Complex& a, const Matrix& b, tree::expression_type t) |
1
|
1697 { |
|
1698 enum RT { RT_unknown, RT_real, RT_complex }; |
|
1699 RT result_type = RT_unknown; |
|
1700 |
|
1701 Matrix result; |
|
1702 ComplexMatrix complex_result; |
143
|
1703 |
1
|
1704 switch (t) |
|
1705 { |
|
1706 case tree::add: |
|
1707 result_type = RT_complex; |
|
1708 complex_result = a + b; |
|
1709 break; |
|
1710 case tree::subtract: |
|
1711 result_type = RT_complex; |
|
1712 complex_result = a - b; |
|
1713 break; |
|
1714 case tree::el_leftdiv: |
|
1715 case tree::leftdiv: |
|
1716 if (a == 0.0) |
|
1717 DIVIDE_BY_ZERO_ERROR; |
164
|
1718 result_type = RT_complex; |
|
1719 complex_result = b / a; |
|
1720 break; |
1
|
1721 case tree::multiply: |
|
1722 case tree::el_mul: |
|
1723 result_type = RT_complex; |
|
1724 complex_result = a * b; |
|
1725 break; |
|
1726 case tree::el_div: |
|
1727 return x_el_div (a, b); |
|
1728 break; |
|
1729 case tree::divide: |
143
|
1730 gripe_nonconformant (1, 1, b.rows (), b.columns ()); |
1
|
1731 break; |
|
1732 case tree::power: |
|
1733 return xpow (a, b); |
|
1734 break; |
|
1735 case tree::elem_pow: |
|
1736 return elem_xpow (a, b); |
|
1737 break; |
|
1738 case tree::cmp_lt: |
|
1739 result_type = RT_real; |
|
1740 result = mx_stupid_bool_op (Matrix_LT, a, b); |
|
1741 break; |
|
1742 case tree::cmp_le: |
|
1743 result_type = RT_real; |
|
1744 result = mx_stupid_bool_op (Matrix_LE, a, b); |
|
1745 break; |
|
1746 case tree::cmp_eq: |
|
1747 result_type = RT_real; |
|
1748 result = mx_stupid_bool_op (Matrix_EQ, a, b); |
|
1749 break; |
|
1750 case tree::cmp_ge: |
|
1751 result_type = RT_real; |
|
1752 result = mx_stupid_bool_op (Matrix_GE, a, b); |
|
1753 break; |
|
1754 case tree::cmp_gt: |
|
1755 result_type = RT_real; |
|
1756 result = mx_stupid_bool_op (Matrix_GT, a, b); |
|
1757 break; |
|
1758 case tree::cmp_ne: |
|
1759 result_type = RT_real; |
|
1760 result = mx_stupid_bool_op (Matrix_NE, a, b); |
|
1761 break; |
|
1762 case tree::and: |
|
1763 result_type = RT_real; |
|
1764 result = mx_stupid_bool_op (Matrix_AND, a, b); |
|
1765 break; |
|
1766 case tree::or: |
|
1767 result_type = RT_real; |
|
1768 result = mx_stupid_bool_op (Matrix_OR, a, b); |
|
1769 break; |
|
1770 default: |
|
1771 panic_impossible (); |
|
1772 break; |
|
1773 } |
|
1774 |
143
|
1775 if (error_state) |
|
1776 return tree_constant (); |
|
1777 |
1
|
1778 assert (result_type != RT_unknown); |
143
|
1779 |
1
|
1780 if (result_type == RT_real) |
|
1781 return tree_constant (result); |
|
1782 else |
|
1783 return tree_constant (complex_result); |
|
1784 } |
|
1785 |
|
1786 /* 11 */ |
|
1787 tree_constant |
164
|
1788 do_binary_op (const Complex& a, const Complex& b, tree::expression_type t) |
1
|
1789 { |
|
1790 enum RT { RT_unknown, RT_real, RT_complex }; |
|
1791 RT result_type = RT_unknown; |
|
1792 |
|
1793 double result = 0.0; |
|
1794 Complex complex_result; |
143
|
1795 |
1
|
1796 switch (t) |
|
1797 { |
|
1798 case tree::add: |
|
1799 result_type = RT_complex; |
|
1800 complex_result = a + b; |
|
1801 break; |
|
1802 case tree::subtract: |
|
1803 result_type = RT_complex; |
|
1804 complex_result = a - b; |
|
1805 break; |
|
1806 case tree::multiply: |
|
1807 case tree::el_mul: |
|
1808 result_type = RT_complex; |
|
1809 complex_result = a * b; |
|
1810 break; |
|
1811 case tree::divide: |
|
1812 case tree::el_div: |
|
1813 result_type = RT_complex; |
|
1814 if (b == 0.0) |
|
1815 DIVIDE_BY_ZERO_ERROR; |
|
1816 complex_result = a / b; |
|
1817 break; |
|
1818 case tree::leftdiv: |
|
1819 case tree::el_leftdiv: |
|
1820 result_type = RT_complex; |
|
1821 if (a == 0.0) |
|
1822 DIVIDE_BY_ZERO_ERROR; |
|
1823 complex_result = b / a; |
|
1824 break; |
|
1825 case tree::power: |
|
1826 case tree::elem_pow: |
|
1827 return xpow (a, b); |
|
1828 break; |
|
1829 case tree::cmp_lt: |
|
1830 result_type = RT_real; |
|
1831 result = real (a) < real (b); |
|
1832 break; |
|
1833 case tree::cmp_le: |
|
1834 result_type = RT_real; |
|
1835 result = real (a) <= real (b); |
|
1836 break; |
|
1837 case tree::cmp_eq: |
|
1838 result_type = RT_real; |
|
1839 result = a == b; |
|
1840 break; |
|
1841 case tree::cmp_ge: |
|
1842 result_type = RT_real; |
|
1843 result = real (a) >= real (b); |
|
1844 break; |
|
1845 case tree::cmp_gt: |
|
1846 result_type = RT_real; |
|
1847 result = real (a) > real (b); |
|
1848 break; |
|
1849 case tree::cmp_ne: |
|
1850 result_type = RT_real; |
|
1851 result = a != b; |
|
1852 break; |
|
1853 case tree::and: |
|
1854 result_type = RT_real; |
|
1855 result = ((a != 0.0) && (b != 0.0)); |
|
1856 break; |
|
1857 case tree::or: |
|
1858 result_type = RT_real; |
|
1859 result = ((a != 0.0) || (b != 0.0)); |
|
1860 break; |
|
1861 default: |
|
1862 panic_impossible (); |
|
1863 break; |
|
1864 } |
|
1865 |
143
|
1866 if (error_state) |
|
1867 return tree_constant (); |
|
1868 |
1
|
1869 assert (result_type != RT_unknown); |
143
|
1870 |
1
|
1871 if (result_type == RT_real) |
|
1872 return tree_constant (result); |
|
1873 else |
|
1874 return tree_constant (complex_result); |
|
1875 } |
|
1876 |
|
1877 /* 12 */ |
|
1878 tree_constant |
164
|
1879 do_binary_op (const Complex& a, const ComplexMatrix& b, |
|
1880 tree::expression_type t) |
1
|
1881 { |
|
1882 enum RT { RT_unknown, RT_real, RT_complex }; |
|
1883 RT result_type = RT_unknown; |
|
1884 |
|
1885 Matrix result; |
|
1886 ComplexMatrix complex_result; |
143
|
1887 |
1
|
1888 switch (t) |
|
1889 { |
|
1890 case tree::add: |
|
1891 result_type = RT_complex; |
|
1892 complex_result = a + b; |
|
1893 break; |
|
1894 case tree::subtract: |
|
1895 result_type = RT_complex; |
|
1896 complex_result = a - b; |
|
1897 break; |
|
1898 case tree::el_leftdiv: |
|
1899 case tree::leftdiv: |
|
1900 if (a == 0.0) |
|
1901 DIVIDE_BY_ZERO_ERROR; |
164
|
1902 result_type = RT_complex; |
|
1903 complex_result = b / a; |
|
1904 break; |
1
|
1905 case tree::multiply: |
|
1906 case tree::el_mul: |
|
1907 result_type = RT_complex; |
|
1908 complex_result = a * b; |
|
1909 break; |
|
1910 case tree::el_div: |
|
1911 return x_el_div (a, b); |
|
1912 break; |
|
1913 case tree::divide: |
143
|
1914 gripe_nonconformant (1, 1, b.rows (), b.columns ()); |
1
|
1915 break; |
|
1916 case tree::power: |
|
1917 return xpow (a, b); |
|
1918 break; |
|
1919 case tree::elem_pow: |
|
1920 return elem_xpow (a, b); |
|
1921 break; |
|
1922 case tree::cmp_lt: |
|
1923 result_type = RT_real; |
|
1924 result = mx_stupid_bool_op (Matrix_LT, a, b); |
|
1925 break; |
|
1926 case tree::cmp_le: |
|
1927 result_type = RT_real; |
|
1928 result = mx_stupid_bool_op (Matrix_LE, a, b); |
|
1929 break; |
|
1930 case tree::cmp_eq: |
|
1931 result_type = RT_real; |
|
1932 result = mx_stupid_bool_op (Matrix_EQ, a, b); |
|
1933 break; |
|
1934 case tree::cmp_ge: |
|
1935 result_type = RT_real; |
|
1936 result = mx_stupid_bool_op (Matrix_GE, a, b); |
|
1937 break; |
|
1938 case tree::cmp_gt: |
|
1939 result_type = RT_real; |
|
1940 result = mx_stupid_bool_op (Matrix_GT, a, b); |
|
1941 break; |
|
1942 case tree::cmp_ne: |
|
1943 result_type = RT_real; |
|
1944 result = mx_stupid_bool_op (Matrix_NE, a, b); |
|
1945 break; |
|
1946 case tree::and: |
|
1947 result_type = RT_real; |
|
1948 result = mx_stupid_bool_op (Matrix_AND, a, b); |
|
1949 break; |
|
1950 case tree::or: |
|
1951 result_type = RT_real; |
|
1952 result = mx_stupid_bool_op (Matrix_OR, a, b); |
|
1953 break; |
|
1954 default: |
|
1955 panic_impossible (); |
|
1956 break; |
|
1957 } |
|
1958 |
143
|
1959 if (error_state) |
|
1960 return tree_constant (); |
|
1961 |
1
|
1962 assert (result_type != RT_unknown); |
143
|
1963 |
1
|
1964 if (result_type == RT_real) |
|
1965 return tree_constant (result); |
|
1966 else |
|
1967 return tree_constant (complex_result); |
|
1968 } |
|
1969 |
|
1970 /* 13 */ |
|
1971 tree_constant |
164
|
1972 do_binary_op (const ComplexMatrix& a, double b, tree::expression_type t) |
1
|
1973 { |
|
1974 enum RT { RT_unknown, RT_real, RT_complex }; |
|
1975 RT result_type = RT_unknown; |
|
1976 |
|
1977 Matrix result; |
|
1978 ComplexMatrix complex_result; |
143
|
1979 |
1
|
1980 switch (t) |
|
1981 { |
|
1982 case tree::add: |
|
1983 result_type = RT_complex; |
|
1984 complex_result = a + b; |
|
1985 break; |
|
1986 case tree::subtract: |
|
1987 result_type = RT_complex; |
|
1988 complex_result = a - b; |
|
1989 break; |
|
1990 case tree::multiply: |
|
1991 case tree::el_mul: |
|
1992 result_type = RT_complex; |
|
1993 complex_result = a * b; |
|
1994 break; |
|
1995 case tree::divide: |
|
1996 case tree::el_div: |
|
1997 result_type = RT_complex; |
|
1998 complex_result = a / b; |
|
1999 break; |
|
2000 case tree::el_leftdiv: |
|
2001 return x_el_div (b, a); |
|
2002 break; |
|
2003 case tree::leftdiv: |
143
|
2004 gripe_nonconformant (a.rows (), a.columns (), 1, 1); |
1
|
2005 break; |
|
2006 case tree::power: |
|
2007 return xpow (a, b); |
|
2008 break; |
|
2009 case tree::elem_pow: |
|
2010 return elem_xpow (a, b); |
|
2011 break; |
|
2012 case tree::cmp_lt: |
|
2013 result_type = RT_real; |
|
2014 result = mx_stupid_bool_op (Matrix_LT, a, b); |
|
2015 break; |
|
2016 case tree::cmp_le: |
|
2017 result_type = RT_real; |
|
2018 result = mx_stupid_bool_op (Matrix_LE, a, b); |
|
2019 break; |
|
2020 case tree::cmp_eq: |
|
2021 result_type = RT_real; |
|
2022 result = mx_stupid_bool_op (Matrix_EQ, a, b); |
|
2023 break; |
|
2024 case tree::cmp_ge: |
|
2025 result_type = RT_real; |
|
2026 result = mx_stupid_bool_op (Matrix_GE, a, b); |
|
2027 break; |
|
2028 case tree::cmp_gt: |
|
2029 result_type = RT_real; |
|
2030 result = mx_stupid_bool_op (Matrix_GT, a, b); |
|
2031 break; |
|
2032 case tree::cmp_ne: |
|
2033 result_type = RT_real; |
|
2034 result = mx_stupid_bool_op (Matrix_NE, a, b); |
|
2035 break; |
|
2036 case tree::and: |
|
2037 result_type = RT_real; |
|
2038 result = mx_stupid_bool_op (Matrix_AND, a, b); |
|
2039 break; |
|
2040 case tree::or: |
|
2041 result_type = RT_real; |
|
2042 result = mx_stupid_bool_op (Matrix_OR, a, b); |
|
2043 break; |
|
2044 default: |
|
2045 panic_impossible (); |
|
2046 break; |
|
2047 } |
|
2048 |
143
|
2049 if (error_state) |
|
2050 return tree_constant (); |
|
2051 |
1
|
2052 assert (result_type != RT_unknown); |
143
|
2053 |
1
|
2054 if (result_type == RT_real) |
|
2055 return tree_constant (result); |
|
2056 else |
|
2057 return tree_constant (complex_result); |
|
2058 } |
|
2059 |
|
2060 /* 14 */ |
|
2061 tree_constant |
164
|
2062 do_binary_op (const ComplexMatrix& a, const Matrix& b, tree::expression_type t) |
1
|
2063 { |
|
2064 enum RT { RT_unknown, RT_real, RT_complex }; |
|
2065 RT result_type = RT_unknown; |
|
2066 |
|
2067 Matrix result; |
|
2068 ComplexMatrix complex_result; |
143
|
2069 |
1
|
2070 switch (t) |
|
2071 { |
|
2072 case tree::add: |
|
2073 result_type = RT_complex; |
|
2074 if (m_add_conform (a, b, 1)) |
|
2075 complex_result = a + b; |
|
2076 break; |
|
2077 case tree::subtract: |
|
2078 result_type = RT_complex; |
|
2079 if (m_add_conform (a, b, 1)) |
|
2080 complex_result = a - b; |
|
2081 break; |
|
2082 case tree::el_mul: |
|
2083 result_type = RT_complex; |
|
2084 if (m_add_conform (a, b, 1)) |
240
|
2085 complex_result = product (a, b); |
1
|
2086 break; |
|
2087 case tree::multiply: |
|
2088 result_type = RT_complex; |
|
2089 if (m_mul_conform (a, b, 1)) |
|
2090 complex_result = a * b; |
|
2091 break; |
|
2092 case tree::el_div: |
|
2093 result_type = RT_complex; |
|
2094 if (m_add_conform (a, b, 1)) |
240
|
2095 complex_result = quotient (a, b); |
1
|
2096 break; |
|
2097 case tree::el_leftdiv: |
|
2098 result_type = RT_complex; |
|
2099 if (m_add_conform (a, b, 1)) |
390
|
2100 complex_result = quotient (b, a); |
1
|
2101 break; |
|
2102 case tree::leftdiv: |
|
2103 return xleftdiv (a, b); |
|
2104 break; |
|
2105 case tree::divide: |
|
2106 return xdiv (a, b); |
|
2107 break; |
|
2108 case tree::power: |
|
2109 error ("can't do A ^ B for A and B both matrices"); |
|
2110 break; |
|
2111 case tree::elem_pow: |
|
2112 if (m_add_conform (a, b, 1)) |
|
2113 return elem_xpow (a, b); |
|
2114 break; |
|
2115 case tree::cmp_lt: |
|
2116 result_type = RT_real; |
|
2117 if (m_add_conform (a, b, 1)) |
|
2118 result = mx_stupid_bool_op (Matrix_LT, a, b); |
|
2119 break; |
|
2120 case tree::cmp_le: |
|
2121 result_type = RT_real; |
|
2122 if (m_add_conform (a, b, 1)) |
|
2123 result = mx_stupid_bool_op (Matrix_LE, a, b); |
|
2124 break; |
|
2125 case tree::cmp_eq: |
|
2126 result_type = RT_real; |
|
2127 if (m_add_conform (a, b, 1)) |
|
2128 result = mx_stupid_bool_op (Matrix_EQ, a, b); |
|
2129 break; |
|
2130 case tree::cmp_ge: |
|
2131 result_type = RT_real; |
|
2132 if (m_add_conform (a, b, 1)) |
|
2133 result = mx_stupid_bool_op (Matrix_GE, a, b); |
|
2134 break; |
|
2135 case tree::cmp_gt: |
|
2136 result_type = RT_real; |
|
2137 if (m_add_conform (a, b, 1)) |
|
2138 result = mx_stupid_bool_op (Matrix_GT, a, b); |
|
2139 break; |
|
2140 case tree::cmp_ne: |
|
2141 result_type = RT_real; |
|
2142 if (m_add_conform (a, b, 1)) |
|
2143 result = mx_stupid_bool_op (Matrix_NE, a, b); |
|
2144 break; |
|
2145 case tree::and: |
|
2146 result_type = RT_real; |
|
2147 if (m_add_conform (a, b, 1)) |
|
2148 result = mx_stupid_bool_op (Matrix_AND, a, b); |
|
2149 break; |
|
2150 case tree::or: |
|
2151 result_type = RT_real; |
|
2152 if (m_add_conform (a, b, 1)) |
|
2153 result = mx_stupid_bool_op (Matrix_OR, a, b); |
|
2154 break; |
|
2155 default: |
|
2156 panic_impossible (); |
|
2157 break; |
|
2158 } |
|
2159 |
143
|
2160 if (error_state) |
|
2161 return tree_constant (); |
|
2162 |
1
|
2163 assert (result_type != RT_unknown); |
143
|
2164 |
1
|
2165 if (result_type == RT_real) |
|
2166 return tree_constant (result); |
|
2167 else |
|
2168 return tree_constant (complex_result); |
|
2169 } |
|
2170 |
|
2171 /* 15 */ |
|
2172 tree_constant |
164
|
2173 do_binary_op (const ComplexMatrix& a, const Complex& b, |
|
2174 tree::expression_type t) |
1
|
2175 { |
|
2176 enum RT { RT_unknown, RT_real, RT_complex }; |
|
2177 RT result_type = RT_unknown; |
|
2178 |
|
2179 Matrix result; |
|
2180 ComplexMatrix complex_result; |
143
|
2181 |
1
|
2182 switch (t) |
|
2183 { |
|
2184 case tree::add: |
|
2185 result_type = RT_complex; |
|
2186 complex_result = a + b; |
|
2187 break; |
|
2188 case tree::subtract: |
|
2189 result_type = RT_complex; |
|
2190 complex_result = a - b; |
|
2191 break; |
|
2192 case tree::multiply: |
|
2193 case tree::el_mul: |
|
2194 result_type = RT_complex; |
|
2195 complex_result = a * b; |
|
2196 break; |
|
2197 case tree::divide: |
|
2198 case tree::el_div: |
|
2199 result_type = RT_complex; |
|
2200 complex_result = a / b; |
|
2201 break; |
|
2202 case tree::el_leftdiv: |
|
2203 return x_el_div (b, a); |
|
2204 break; |
|
2205 case tree::leftdiv: |
143
|
2206 gripe_nonconformant (a.rows (), a.columns (), 1, 1); |
1
|
2207 break; |
|
2208 case tree::power: |
|
2209 return xpow (a, b); |
|
2210 break; |
|
2211 case tree::elem_pow: |
|
2212 return elem_xpow (a, b); |
|
2213 break; |
|
2214 case tree::cmp_lt: |
|
2215 result_type = RT_real; |
|
2216 result = mx_stupid_bool_op (Matrix_LT, a, b); |
|
2217 break; |
|
2218 case tree::cmp_le: |
|
2219 result_type = RT_real; |
|
2220 result = mx_stupid_bool_op (Matrix_LE, a, b); |
|
2221 break; |
|
2222 case tree::cmp_eq: |
|
2223 result_type = RT_real; |
|
2224 result = mx_stupid_bool_op (Matrix_EQ, a, b); |
|
2225 break; |
|
2226 case tree::cmp_ge: |
|
2227 result_type = RT_real; |
|
2228 result = mx_stupid_bool_op (Matrix_GE, a, b); |
|
2229 break; |
|
2230 case tree::cmp_gt: |
|
2231 result_type = RT_real; |
|
2232 result = mx_stupid_bool_op (Matrix_GT, a, b); |
|
2233 break; |
|
2234 case tree::cmp_ne: |
|
2235 result_type = RT_real; |
|
2236 result = mx_stupid_bool_op (Matrix_NE, a, b); |
|
2237 break; |
|
2238 case tree::and: |
|
2239 result_type = RT_real; |
|
2240 result = mx_stupid_bool_op (Matrix_AND, a, b); |
|
2241 break; |
|
2242 case tree::or: |
|
2243 result_type = RT_real; |
|
2244 result = mx_stupid_bool_op (Matrix_OR, a, b); |
|
2245 break; |
|
2246 default: |
|
2247 panic_impossible (); |
|
2248 break; |
|
2249 } |
|
2250 |
143
|
2251 if (error_state) |
|
2252 return tree_constant (); |
|
2253 |
1
|
2254 assert (result_type != RT_unknown); |
143
|
2255 |
1
|
2256 if (result_type == RT_real) |
|
2257 return tree_constant (result); |
|
2258 else |
|
2259 return tree_constant (complex_result); |
|
2260 } |
|
2261 |
|
2262 /* 16 */ |
|
2263 tree_constant |
164
|
2264 do_binary_op (const ComplexMatrix& a, const ComplexMatrix& b, |
|
2265 tree::expression_type t) |
1
|
2266 { |
|
2267 enum RT { RT_unknown, RT_real, RT_complex }; |
|
2268 RT result_type = RT_unknown; |
|
2269 |
|
2270 Matrix result; |
|
2271 ComplexMatrix complex_result; |
143
|
2272 |
1
|
2273 switch (t) |
|
2274 { |
|
2275 case tree::add: |
|
2276 result_type = RT_complex; |
|
2277 if (m_add_conform (a, b, 1)) |
|
2278 complex_result = a + b; |
|
2279 break; |
|
2280 case tree::subtract: |
|
2281 result_type = RT_complex; |
|
2282 if (m_add_conform (a, b, 1)) |
|
2283 complex_result = a - b; |
|
2284 break; |
|
2285 case tree::el_mul: |
|
2286 result_type = RT_complex; |
|
2287 if (m_add_conform (a, b, 1)) |
240
|
2288 complex_result = product (a, b); |
1
|
2289 break; |
|
2290 case tree::multiply: |
|
2291 result_type = RT_complex; |
|
2292 if (m_mul_conform (a, b, 1)) |
|
2293 complex_result = a * b; |
|
2294 break; |
|
2295 case tree::el_div: |
|
2296 result_type = RT_complex; |
|
2297 if (m_add_conform (a, b, 1)) |
240
|
2298 complex_result = quotient (a, b); |
1
|
2299 break; |
|
2300 case tree::el_leftdiv: |
|
2301 result_type = RT_complex; |
|
2302 if (m_add_conform (a, b, 1)) |
240
|
2303 complex_result = quotient (b, a); |
1
|
2304 break; |
|
2305 case tree::leftdiv: |
|
2306 return xleftdiv (a, b); |
|
2307 break; |
|
2308 case tree::divide: |
|
2309 return xdiv (a, b); |
|
2310 break; |
|
2311 case tree::power: |
|
2312 error ("can't do A ^ B for A and B both matrices"); |
|
2313 break; |
|
2314 case tree::elem_pow: |
|
2315 if (m_add_conform (a, b, 1)) |
|
2316 return elem_xpow (a, b); |
|
2317 break; |
|
2318 case tree::cmp_lt: |
|
2319 result_type = RT_real; |
|
2320 if (m_add_conform (a, b, 1)) |
|
2321 result = mx_stupid_bool_op (Matrix_LT, a, b); |
|
2322 break; |
|
2323 case tree::cmp_le: |
|
2324 result_type = RT_real; |
|
2325 if (m_add_conform (a, b, 1)) |
|
2326 result = mx_stupid_bool_op (Matrix_LE, a, b); |
|
2327 break; |
|
2328 case tree::cmp_eq: |
|
2329 result_type = RT_real; |
|
2330 if (m_add_conform (a, b, 1)) |
|
2331 result = mx_stupid_bool_op (Matrix_EQ, a, b); |
|
2332 break; |
|
2333 case tree::cmp_ge: |
|
2334 result_type = RT_real; |
|
2335 if (m_add_conform (a, b, 1)) |
|
2336 result = mx_stupid_bool_op (Matrix_GE, a, b); |
|
2337 break; |
|
2338 case tree::cmp_gt: |
|
2339 result_type = RT_real; |
|
2340 if (m_add_conform (a, b, 1)) |
|
2341 result = mx_stupid_bool_op (Matrix_GT, a, b); |
|
2342 break; |
|
2343 case tree::cmp_ne: |
|
2344 result_type = RT_real; |
|
2345 if (m_add_conform (a, b, 1)) |
|
2346 result = mx_stupid_bool_op (Matrix_NE, a, b); |
|
2347 break; |
|
2348 case tree::and: |
|
2349 result_type = RT_real; |
|
2350 if (m_add_conform (a, b, 1)) |
|
2351 result = mx_stupid_bool_op (Matrix_AND, a, b); |
|
2352 break; |
|
2353 case tree::or: |
|
2354 result_type = RT_real; |
|
2355 if (m_add_conform (a, b, 1)) |
|
2356 result = mx_stupid_bool_op (Matrix_OR, a, b); |
|
2357 break; |
|
2358 default: |
|
2359 panic_impossible (); |
|
2360 break; |
|
2361 } |
|
2362 |
143
|
2363 if (error_state) |
|
2364 return tree_constant (); |
|
2365 |
1
|
2366 assert (result_type != RT_unknown); |
143
|
2367 |
1
|
2368 if (result_type == RT_real) |
|
2369 return tree_constant (result); |
|
2370 else |
|
2371 return tree_constant (complex_result); |
|
2372 } |
|
2373 |
|
2374 /* |
|
2375 ;;; Local Variables: *** |
|
2376 ;;; mode: C++ *** |
|
2377 ;;; page-delimiter: "^/\\*" *** |
|
2378 ;;; End: *** |
|
2379 */ |