X-Git-Url: http://nsz.repo.hu/git/?a=blobdiff_plain;f=ir%2Fadt%2Fgaussjordan.c;h=34cd9b64669345f838bfaa00c61db9dc5964b925;hb=ebdaa596d904a0651c8d794481288d2d629bdb3a;hp=678e43b948fbbf59748b3b7a5c17e9ebeeb8921a;hpb=b78bdd4d94de46de4156272e6dbfe44e97933a5b;p=libfirm diff --git a/ir/adt/gaussjordan.c b/ir/adt/gaussjordan.c index 678e43b94..34cd9b646 100644 --- a/ir/adt/gaussjordan.c +++ b/ir/adt/gaussjordan.c @@ -32,10 +32,7 @@ /* returns 0 if successful */ /* returns -1 if ill-conditioned matrix */ /*------------------------------------------------------*/ - -#ifdef HAVE_CONFIG_H #include "config.h" -#endif #include #include @@ -45,91 +42,96 @@ int firm_gaussjordansolve(double *A, double *vec, int nsize) { - int i, j, row, col, col2, biggest_r, biggest_c; - double big, temp, sum; - double *x = xmalloc(nsize * sizeof(double)); - double *scramvec = xmalloc(nsize * sizeof(double)); + int i, j, row, col, col2, biggest_r = 0, biggest_c = 0, t; + double big, temp, sum; + double *scramvec = XMALLOCN(double, nsize); + int *x = XMALLOCN(int, nsize); + int res = 0; #define _A(row,col) A[(row)*nsize + (col)] - /* init x[] */ - for(i=0;i big) { - biggest_r = row; - biggest_c = col2; - big = temp; /* largest element left */ + /* init x[] */ + for (i = 0; i < nsize; ++i) + x[i] = i; + + /* triangularize A */ + /* ie A has zeros below it's diagonal */ + for (col = 0; col < nsize - 1; ++col) { + big = 0; + /* find the largest left in LRH box */ + for (row = col; row < nsize; ++row) { + for (col2 = col; col2 < nsize; ++col2) { + temp = fabs(_A(row,col2)); + if (temp > big) { + biggest_r = row; + biggest_c = col2; + big = temp; /* largest element left */ + } + } + } + if (big < SMALL) { + res = -1; + goto end; + } + + /* swap rows */ + for(i=0;i=0;i--) { + sum = 0; + for(j=i+1;j=0;i--) { - sum = 0; - for(j=i+1;j