1 /* GENIUS Calculator
2 * Copyright (C) 1997-2012 Jiri (George) Lebl
3 *
4 * Author: Jiri (George) Lebl
5 *
6 * This file is part of Genius.
7 *
8 * Genius is free software: you can redistribute it and/or modify
9 * it under the terms of the GNU General Public License as published by
10 * the Free Software Foundation, either version 3 of the License, or
11 * (at your option) any later version.
12 *
13 * This program is distributed in the hope that it will be useful,
14 * but WITHOUT ANY WARRANTY; without even the implied warranty of
15 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
16 * GNU General Public License for more details.
17 *
18 * You should have received a copy of the GNU General Public License
19 * along with this program. If not, see <http://www.gnu.org/licenses/>.
20 */
21
22 #include "config.h"
23
24 #include <stdio.h>
25 #include <string.h>
26 #include <glib.h>
27 #include "calc.h"
28 #include "mpwrap.h"
29 #include "eval.h"
30 #include "dict.h"
31 #include "util.h"
32 #include "funclib.h"
33 #include "matrix.h"
34 #include "matrixw.h"
35
36 #include "matop.h"
37
38 gboolean
gel_is_matrix_value_only(GelMatrixW * m)39 gel_is_matrix_value_only (GelMatrixW *m)
40 {
41 int i, j, w, h;
42 if (m->cached_value_only)
43 return m->value_only;
44 w = gel_matrixw_width (m);
45 h = gel_matrixw_height (m);
46 for (j = 0; j < h; j++) {
47 for (i = 0; i < w; i++) {
48 GelETree *n = gel_matrixw_get_index(m,i,j);
49 if(n && n->type != GEL_VALUE_NODE) {
50 m->cached_value_only = 1;
51 m->value_only = 0;
52 m->cached_value_only_real = 1;
53 m->value_only_real = 0;
54 m->cached_value_only_rational = 1;
55 m->value_only_rational = 0;
56 m->cached_value_only_integer = 1;
57 m->value_only_integer = 0;
58 return FALSE;
59 }
60 }
61 }
62 m->cached_value_only = 1;
63 m->value_only = 1;
64 return TRUE;
65 }
66
67 gboolean
gel_is_matrix_value_or_bool_only(GelMatrixW * m)68 gel_is_matrix_value_or_bool_only (GelMatrixW *m)
69 {
70 int i, j, w, h;
71 gboolean got_bools = FALSE;
72 if (m->cached_value_or_bool_only)
73 return m->value_or_bool_only;
74 w = gel_matrixw_width (m);
75 h = gel_matrixw_height (m);
76 for (j = 0; j < h; j++) {
77 for (i = 0; i < w; i++) {
78 GelETree *n = gel_matrixw_get_index(m,i,j);
79 if ( ! n)
80 continue;
81 if (n->type == GEL_BOOL_NODE) {
82 got_bools = TRUE;
83 continue;
84 }
85
86 if (n->type != GEL_VALUE_NODE) {
87 m->cached_value_or_bool_only = 1;
88 m->value_or_bool_only = 0;
89 m->cached_value_only = 1;
90 m->value_only = 0;
91 m->cached_value_only_real = 1;
92 m->value_only_real = 0;
93 m->cached_value_only_rational = 1;
94 m->value_only_rational = 0;
95 m->cached_value_only_integer = 1;
96 m->value_only_integer = 0;
97 return FALSE;
98 }
99 }
100 }
101 m->cached_value_or_bool_only = 1;
102 m->value_or_bool_only = 1;
103
104 m->cached_value_only = 1;
105 if (got_bools)
106 m->value_only = 0;
107 else
108 m->value_only = 1;
109 return TRUE;
110 }
111
112
113 gboolean
gel_is_matrix_value_only_real(GelMatrixW * m)114 gel_is_matrix_value_only_real (GelMatrixW *m)
115 {
116 int i, j, w, h;
117 if (m->cached_value_only_real)
118 return m->value_only_real;
119 w = gel_matrixw_width (m);
120 h = gel_matrixw_height (m);
121 for (j = 0; j < h; j++) {
122 for (i = 0; i < w; i++) {
123 GelETree *n = gel_matrixw_get_index(m,i,j);
124 if (n != NULL &&
125 (n->type != GEL_VALUE_NODE ||
126 mpw_is_complex (n->val.value))) {
127 m->cached_value_only_real = 1;
128 m->value_only_real = 0;
129 return FALSE;
130 }
131 }
132 }
133 m->cached_value_only = 1;
134 m->value_only = 1;
135 m->cached_value_only_real = 1;
136 m->value_only_real = 1;
137 return TRUE;
138 }
139
140 gboolean
gel_is_matrix_value_only_rational(GelMatrixW * m)141 gel_is_matrix_value_only_rational (GelMatrixW *m)
142 {
143 int i, j, w, h;
144 if (m->cached_value_only_rational)
145 return m->value_only_rational;
146 w = gel_matrixw_width (m);
147 h = gel_matrixw_height (m);
148 for (j = 0; j < h; j++) {
149 for (i = 0; i < w; i++) {
150 GelETree *n = gel_matrixw_get_index(m,i,j);
151 if (n != NULL &&
152 (n->type != GEL_VALUE_NODE ||
153 mpw_is_complex (n->val.value) ||
154 mpw_is_real_part_float (n->val.value))) {
155 m->cached_value_only_rational = 1;
156 m->value_only_rational = 0;
157 return FALSE;
158 }
159 }
160 }
161 m->cached_value_only = 1;
162 m->value_only = 1;
163 m->cached_value_only_rational = 1;
164 m->value_only_rational = 1;
165 return TRUE;
166 }
167
168 gboolean
gel_is_matrix_value_only_integer(GelMatrixW * m)169 gel_is_matrix_value_only_integer (GelMatrixW *m)
170 {
171 int i, j, w, h;
172 if (m->cached_value_only_integer)
173 return m->value_only_integer;
174 w = gel_matrixw_width (m);
175 h = gel_matrixw_height (m);
176 for (j = 0; j < h; j++) {
177 for (i = 0; i < w; i++) {
178 GelETree *n = gel_matrixw_get_index(m,i,j);
179 if (n != NULL &&
180 (n->type != GEL_VALUE_NODE ||
181 mpw_is_complex (n->val.value) ||
182 ! mpw_is_integer (n->val.value))) {
183 m->cached_value_only_integer = 1;
184 m->value_only_integer = 0;
185 return FALSE;
186 }
187 }
188 }
189 m->cached_value_only = 1;
190 m->value_only = 1;
191 m->cached_value_only_rational = 1;
192 m->value_only_rational = 1;
193 m->cached_value_only_integer = 1;
194 m->value_only_integer = 1;
195 return TRUE;
196 }
197
198 void
gel_matrix_conjugate_transpose(GelMatrixW * m)199 gel_matrix_conjugate_transpose (GelMatrixW *m)
200 {
201 int i, j, w, h;
202
203 if (gel_is_matrix_value_only_real (m)) {
204 m->tr = !m->tr;
205 return;
206 }
207
208 gel_matrixw_make_private (m, FALSE /* kill_type_caches */);
209 m->tr = !m->tr;
210 w = gel_matrixw_width (m);
211 h = gel_matrixw_height (m);
212 for (j = 0; j < h; j++) {
213 for (i = 0; i < w; i++) {
214 GelETree *n = gel_matrixw_get_index (m, i, j);
215 if (n == NULL)
216 continue;
217 if (n->type == GEL_VALUE_NODE) {
218 if (mpw_is_complex (n->val.value))
219 mpw_conj (n->val.value, n->val.value);
220 } else {
221 GelETree *nn;
222 GEL_GET_NEW_NODE (nn);
223 nn->type = GEL_OPERATOR_NODE;
224 nn->op.oper = GEL_E_DIRECTCALL;
225 nn->op.nargs = 2;
226
227 GEL_GET_NEW_NODE (nn->op.args);
228 nn->op.args->type = GEL_IDENTIFIER_NODE;
229 nn->op.args->id.id = d_intern ("conj");
230 nn->op.args->id.uninitialized = FALSE;
231
232 nn->op.args->any.next = n;
233 n->any.next = NULL;
234 gel_matrixw_set_index (m, i, j) = nn;
235 }
236 }
237 }
238 }
239
240 void
gel_value_matrix_multiply(GelMatrixW * res,GelMatrixW * m1,GelMatrixW * m2,mpw_ptr modulo)241 gel_value_matrix_multiply (GelMatrixW *res, GelMatrixW *m1, GelMatrixW *m2,
242 mpw_ptr modulo)
243 {
244 int i, j, k, w, h, m1w;
245 mpw_t tmp;
246 mpw_init(tmp);
247 gel_matrixw_make_private(res, TRUE /* kill_type_caches */);
248
249 w = gel_matrixw_width (res);
250 h = gel_matrixw_height (res);
251 m1w = gel_matrixw_width (m1);
252 for (j = 0; j < h; j++) {
253 for (i = 0; i < w; i++) {
254 gboolean got_something = FALSE;
255 mpw_t accu;
256 mpw_init(accu);
257 for (k = 0; k < m1w; k++) {
258 GelETree *l = gel_matrixw_get_index(m1,k,j);
259 GelETree *r = gel_matrixw_get_index(m2,i,k);
260
261 /* if both zero add nothing */
262 if (l == NULL || r == NULL)
263 continue;
264
265 got_something = TRUE;
266
267 mpw_mul(tmp,l->val.value,r->val.value);
268 mpw_add(accu,accu,tmp);
269 if (modulo != NULL) {
270 mpw_mod (accu, accu, modulo);
271 if G_UNLIKELY (gel_error_num != 0) { /*FIXME: for now ignore errors in moding*/
272 gel_error_num = 0;
273 }
274 if (mpw_sgn (accu) < 0)
275 mpw_add (accu, modulo, accu);
276 }
277 /*XXX: are there any problems that could occur
278 here? ... I don't seem to see any, if there
279 are catch them here*/
280 }
281 if (got_something) {
282 gel_matrixw_set_index(res,i,j) = gel_makenum_use(accu);
283 } else {
284 gel_matrixw_set_index(res,i,j) = NULL;
285 mpw_clear (accu);
286 }
287 }
288 }
289 mpw_clear(tmp);
290 }
291
292 /* m must be made private before */
293 static void
swap_rows(GelMatrixW * m,int row,int row2)294 swap_rows(GelMatrixW *m, int row, int row2)
295 {
296 int i, w;
297 if(row==row2) return;
298
299 /* no need for this, only used within gauss and matrix is already private
300 gel_matrixw_make_private(m);*/
301
302 w = gel_matrixw_width (m);
303 for (i = 0; i < w; i++) {
304 GelETree *t = gel_matrixw_get_index(m,i,row);
305 gel_matrixw_set_index(m,i,row) = gel_matrixw_get_index(m,i,row2);
306 gel_matrixw_set_index(m,i,row2) = t;
307 }
308 }
309
310 /* m must be made private before */
311 static gboolean
div_row(GelCtx * ctx,GelMatrixW * m,int row,mpw_t div)312 div_row (GelCtx *ctx, GelMatrixW *m, int row, mpw_t div)
313 {
314 int i, w;
315 gboolean ret = TRUE;
316
317 if (mpw_eql_ui (div, 1))
318 return TRUE;
319
320 /* no need for this, only used within gauss and matrix is already private
321 gel_matrixw_make_private(m);*/
322
323 w = gel_matrixw_width (m);
324 for (i = 0; i < w; i++) {
325 GelETree *t = gel_matrixw_get_index(m,i,row);
326 if(t) {
327 mpw_div(t->val.value,t->val.value,div);
328 if (ctx->modulo != NULL) {
329 gel_mod_node (ctx, t);
330 /* can't mod so we have a singular matrix / system */
331 if (t != NULL && t->type != GEL_VALUE_NODE)
332 ret = FALSE;
333 }
334 }
335 }
336 return ret;
337 }
338
339 /* m must be made private before */
340 static gboolean
mul_sub_row(GelCtx * ctx,GelMatrixW * m,int row,mpw_t mul,int row2)341 mul_sub_row (GelCtx *ctx, GelMatrixW *m, int row, mpw_t mul, int row2)
342 {
343 int i, w;
344 static mpw_t tmp;
345 static gboolean tmp_inited = FALSE;
346 gboolean ret = TRUE;
347
348 /* no need for this, only used within gauss and matrix is already private
349 gel_matrixw_make_private(m);*/
350
351 if G_UNLIKELY ( ! tmp_inited) {
352 mpw_init(tmp);
353 tmp_inited = TRUE;
354 }
355 w = gel_matrixw_width(m);
356 for (i = 0; i < w; i++) {
357 GelETree *t = gel_matrixw_get_index(m,i,row);
358 if (t && ! mpw_zero_p (t->val.value)) {
359 GelETree *t2 = gel_matrixw_get_index(m,i,row2);
360 mpw_mul(tmp,t->val.value,mul);
361 if (t2 == NULL) {
362 mpw_neg(tmp,tmp);
363 t2 = gel_makenum_use(tmp);
364 gel_matrixw_set_index(m,i,row2) = t2;
365 mpw_init(tmp);
366 } else if ( ! mpw_is_complex_float (tmp) &&
367 mpw_symbolic_eql (t2->val.value, tmp)) {
368 gel_freetree (t2);
369 gel_matrixw_set_index(m,i,row2) = NULL;
370 } else {
371 mpw_sub (t2->val.value,
372 t2->val.value, tmp);
373 }
374 if (ctx->modulo != NULL && t2 != NULL) {
375 gel_mod_node (ctx, t2);
376 /* can't mod so we have a singular matrix / system */
377 if (t2 != NULL && t2->type != GEL_VALUE_NODE)
378 ret = FALSE;
379 }
380 }
381 }
382 return ret;
383 }
384
385 /*NOTE: if simul is passed then we assume that it's the same size as m*/
386 /* return FALSE if singular */
387 /* FIXME: if modular arithmetic is on, work over the modulo properly!!!! */
388 gboolean
gel_value_matrix_gauss(GelCtx * ctx,GelMatrixW * m,gboolean reduce,gboolean uppertriang,gboolean stopsing,gboolean stopnonsing,mpw_ptr detop,GelMatrixW * simul)389 gel_value_matrix_gauss (GelCtx *ctx,
390 GelMatrixW *m,
391 gboolean reduce,
392 gboolean uppertriang,
393 gboolean stopsing,
394 gboolean stopnonsing,
395 mpw_ptr detop,
396 GelMatrixW *simul)
397 {
398 int i, j, d, ii, w, h;
399 GelETree *piv;
400 mpw_t tmp;
401 gboolean ret = TRUE;
402 gboolean made_private = FALSE;
403 gboolean matrix_rational = FALSE;
404 int *pivots = NULL;
405 int pivots_max = -1;
406
407 w = gel_matrixw_width (m);
408 h = gel_matrixw_height (m);
409
410 if(detop)
411 mpw_set_ui(detop,1);
412
413 if (m->rref) {
414 /* test for singularity */
415 if (w > h) {
416 ret = FALSE;
417 } else {
418 GelETree *t = gel_matrixw_get_indexii(m,w-1);
419 if (t == NULL ||
420 mpw_zero_p (t->val.value))
421 ret = FALSE;
422 }
423 return ret;
424 }
425
426 matrix_rational = gel_is_matrix_value_only_rational (m);
427
428 mpw_init(tmp);
429 d = 0;
430
431 if (reduce) {
432 pivots = g_alloca (sizeof(int) * h);
433 }
434
435 for (i = 0; i < w && d < h; i++) {
436 if (matrix_rational) {
437 for (j = d; j < h; j++) {
438 GelETree *t = gel_matrixw_get_index(m,i,j);
439 if (t != NULL &&
440 ! mpw_zero_p (t->val.value))
441 break;
442 }
443 } else {
444 /* kind of a hack */
445 int bestj = h;
446 mpw_t best_abs_sq;
447 GelETree *bestpiv = NULL;
448
449 mpw_init (best_abs_sq);
450 for (j = d; j < h; j++) {
451 GelETree *t = gel_matrixw_get_index(m,i,j);
452 if (t != NULL &&
453 ! mpw_zero_p (t->val.value)) {
454 if (bestpiv == NULL) {
455 bestpiv = t;
456 bestj = j;
457 } else {
458 mpw_abs_sq (tmp, t->val.value);
459 if (mpw_cmp (tmp, best_abs_sq) > 0) {
460 bestpiv = t;
461 bestj = j;
462 }
463 }
464 }
465 }
466 mpw_clear (best_abs_sq);
467
468 j = bestj;
469 }
470
471 if (j == h) {
472 if(stopsing) {
473 mpw_clear(tmp);
474 return FALSE;
475 }
476 continue;
477 }
478
479 if ( ! made_private) {
480 gel_matrixw_make_private (m, TRUE /* kill_type_caches */);
481 if (simul)
482 gel_matrixw_make_private (simul, TRUE /* kill_type_caches */);
483 made_private = TRUE;
484
485 /* the matrix will be value only */
486 m->cached_value_only = 1;
487 m->value_only = 1;
488
489 if (matrix_rational) {
490 /* the matrix will be value only rational */
491 m->cached_value_only_rational = 1;
492 m->value_only_rational = 1;
493 }
494 }
495
496 if (j > d) {
497 swap_rows(m,j,d);
498 if(simul)
499 swap_rows(simul,j,d);
500 if(detop)
501 mpw_neg(detop,detop);
502 }
503
504 piv = gel_matrixw_index(m,i,d);
505
506 for (j = d+1; j < h; j++) {
507 GelETree *t = gel_matrixw_get_index(m,i,j);
508 /* Assume t is already reduced mod ctx->modulo
509 * if appropriate */
510 if (t != NULL &&
511 ! mpw_zero_p (t->val.value)) {
512 mpw_div(tmp,t->val.value,piv->val.value);
513 if ( ! mul_sub_row (ctx, m, d, tmp, j) &&
514 stopsing) {
515 mpw_clear(tmp);
516 return FALSE;
517 }
518 if(simul) {
519 if ( ! mul_sub_row (ctx, simul, d, tmp, j) &&
520 stopsing) {
521 mpw_clear(tmp);
522 return FALSE;
523 }
524 }
525 }
526 }
527
528
529 /*we just want to do an upper triangular matrix*/
530 if(uppertriang) {
531 d++;
532 continue;
533 }
534
535 if(reduce) {
536 pivots[d] = i;
537 pivots_max = d;
538 }
539
540 /* make pivot 1 */
541 for (ii = i+1; ii < w; ii++) {
542 GelETree *t = gel_matrixw_get_index(m,ii,d);
543 if(t) {
544 mpw_div(t->val.value,
545 t->val.value,
546 piv->val.value);
547 if (ctx->modulo != NULL) {
548 gel_mod_node (ctx, t);
549 if (stopsing &&
550 t != NULL &&
551 t->type != GEL_VALUE_NODE) {
552 mpw_clear(tmp);
553 return FALSE;
554 }
555 }
556 }
557 }
558 if(detop)
559 mpw_div (detop, detop, piv->val.value);
560 if(simul) {
561 if ( ! div_row (ctx, simul, d, piv->val.value) &&
562 stopsing) {
563 mpw_clear(tmp);
564 return FALSE;
565 }
566 }
567 mpw_set_ui(piv->val.value,1);
568 d++;
569 }
570
571 if (d < w)
572 ret = FALSE;
573
574 if (stopnonsing && d == w) {
575 mpw_clear(tmp);
576 return TRUE;
577 }
578
579 if(reduce) {
580 for(d = pivots_max; d >= 0; d--) {
581 i = pivots[d];
582 for(j=0;j<d;j++) {
583 GelETree *t = gel_matrixw_get_index(m,i,j);
584 if (t != NULL &&
585 ! mpw_zero_p (t->val.value)) {
586 /* subtle: must copy t->val.value,
587 * else we wipe it out */
588 mpw_set (tmp, t->val.value);
589 if ( ! mul_sub_row (ctx, m, d, tmp, j) &&
590 stopsing) {
591 mpw_clear(tmp);
592 return FALSE;
593 }
594 if(simul) {
595 if ( ! mul_sub_row (ctx, simul, d, tmp, j) &&
596 stopsing) {
597 mpw_clear(tmp);
598 return FALSE;
599 }
600 }
601 }
602 }
603 }
604 }
605
606 if (detop && ctx->modulo != NULL) {
607 /* FIXME: this may fail!!! */
608 gel_mod_integer_rational (detop, ctx->modulo);
609 }
610
611 if (reduce && ! uppertriang)
612 m->rref = 1;
613
614 mpw_clear(tmp);
615 return ret;
616 }
617
618
619 static void
det2(mpw_t rop,GelMatrixW * m)620 det2(mpw_t rop, GelMatrixW *m)
621 {
622 mpw_t tmp;
623 mpw_init(tmp);
624 mpw_mul(rop,gel_matrixw_index(m,0,0)->val.value,
625 gel_matrixw_index(m,1,1)->val.value);
626 mpw_mul(tmp,gel_matrixw_index(m,1,0)->val.value,
627 gel_matrixw_index(m,0,1)->val.value);
628 mpw_sub(rop,rop,tmp);
629 mpw_clear(tmp);
630 }
631
632 static void
det3(mpw_t rop,GelMatrixW * m)633 det3(mpw_t rop, GelMatrixW *m)
634 {
635 mpw_t tmp;
636 mpw_init(tmp);
637
638 mpw_mul(rop,gel_matrixw_index(m,0,0)->val.value,
639 gel_matrixw_index(m,1,1)->val.value);
640 mpw_mul(rop,rop,
641 gel_matrixw_index(m,2,2)->val.value);
642
643 mpw_mul(tmp,gel_matrixw_index(m,1,0)->val.value,
644 gel_matrixw_index(m,2,1)->val.value);
645 mpw_mul(tmp,tmp,
646 gel_matrixw_index(m,0,2)->val.value);
647 mpw_add(rop,rop,tmp);
648
649 mpw_mul(tmp,gel_matrixw_index(m,2,0)->val.value,
650 gel_matrixw_index(m,0,1)->val.value);
651 mpw_mul(tmp,tmp,
652 gel_matrixw_index(m,1,2)->val.value);
653 mpw_add(rop,rop,tmp);
654
655 mpw_mul(tmp,gel_matrixw_index(m,2,0)->val.value,
656 gel_matrixw_index(m,1,1)->val.value);
657 mpw_mul(tmp,tmp,
658 gel_matrixw_index(m,0,2)->val.value);
659 mpw_sub(rop,rop,tmp);
660
661 mpw_mul(tmp,gel_matrixw_index(m,1,0)->val.value,
662 gel_matrixw_index(m,0,1)->val.value);
663 mpw_mul(tmp,tmp,
664 gel_matrixw_index(m,2,2)->val.value);
665 mpw_sub(rop,rop,tmp);
666
667 mpw_mul(tmp,gel_matrixw_index(m,0,0)->val.value,
668 gel_matrixw_index(m,2,1)->val.value);
669 mpw_mul(tmp,tmp,
670 gel_matrixw_index(m,1,2)->val.value);
671 mpw_sub(rop,rop,tmp);
672
673 mpw_clear(tmp);
674 }
675
676 gboolean
gel_value_matrix_det(GelCtx * ctx,mpw_t rop,GelMatrixW * m)677 gel_value_matrix_det (GelCtx *ctx, mpw_t rop, GelMatrixW *m)
678 {
679 int w = gel_matrixw_width(m);
680 int h = gel_matrixw_height(m);
681 GelMatrixW *mm;
682 mpw_t tmp;
683 int i;
684
685 /* only works properly without modulo, but modulo is taken
686 * care of after the det function is executed */
687 g_assert (ctx->modulo == NULL);
688
689 if(w != h) {
690 gel_errorout (_("Determinant of a non-square matrix is undefined"));
691 return FALSE;
692 }
693
694 /* If we already are in rref form just compute determinant */
695 if (m->rref) {
696 mpw_set_ui (rop, 1);
697 for (i = 0; i < w; i++) {
698 GelETree *t = gel_matrixw_get_indexii (m, i);
699 if (t == NULL ||
700 mpw_zero_p (t->val.value)) {
701 mpw_set_ui (rop, 0);
702 return TRUE;
703 }
704 /* row reduced form has 1's on the diagonal! */
705 /*mpw_mul(rop,rop,t->val.value);*/
706 }
707 return TRUE;
708 }
709
710
711 switch(w) {
712 case 1:
713 mpw_set(rop,gel_matrixw_index(m,0,0)->val.value);
714 break;
715 case 2:
716 det2(rop,m);
717 break;
718 case 3:
719 det3(rop,m);
720 break;
721 default:
722 mpw_init(tmp);
723 mm = gel_matrixw_copy(m);
724 gel_value_matrix_gauss(ctx,mm,
725 FALSE /* reduce */,
726 TRUE /* uppertriang */,
727 FALSE /* stopsing */,
728 FALSE /* stopnonsing */,
729 tmp,
730 NULL);
731 mpw_mul(rop,tmp,gel_matrixw_index(mm,0,0)->val.value);
732 mpw_clear(tmp);
733 for (i = 1; i < w; i++) {
734 GelETree *t = gel_matrixw_get_indexii(mm,i);
735 if (t == NULL) {
736 gel_matrixw_free(mm);
737 mpw_set_ui(rop,0);
738 return TRUE;
739 }
740 mpw_mul(rop,rop,t->val.value);
741 }
742 gel_matrixw_free(mm);
743 break;
744 }
745 return TRUE;
746 }
747