[Scipy-svn] r6355 - in trunk/scipy/sparse/linalg/dsolve: . SuperLU SuperLU/SRC

scipy-svn at scipy.org scipy-svn at scipy.org
Tue Apr 27 17:57:39 EDT 2010


Author: ptvirtan
Date: 2010-04-27 16:57:39 -0500 (Tue, 27 Apr 2010)
New Revision: 6355

Modified:
   trunk/scipy/sparse/linalg/dsolve/SuperLU/SRC/ilu_cpivotL.c
   trunk/scipy/sparse/linalg/dsolve/SuperLU/SRC/ilu_dpivotL.c
   trunk/scipy/sparse/linalg/dsolve/SuperLU/SRC/ilu_spivotL.c
   trunk/scipy/sparse/linalg/dsolve/SuperLU/SRC/ilu_zpivotL.c
   trunk/scipy/sparse/linalg/dsolve/SuperLU/scipychanges.txt
   trunk/scipy/sparse/linalg/dsolve/_superluobject.c
Log:
BUG: sparse.linalg.dsolve/SuperLU: patch SuperLU sources to not exit(1) when ILU decomposition encounters singularity; instead, raise a Python exception

Modified: trunk/scipy/sparse/linalg/dsolve/SuperLU/SRC/ilu_cpivotL.c
===================================================================
--- trunk/scipy/sparse/linalg/dsolve/SuperLU/SRC/ilu_cpivotL.c	2010-04-27 21:57:15 UTC (rev 6354)
+++ trunk/scipy/sparse/linalg/dsolve/SuperLU/SRC/ilu_cpivotL.c	2010-04-27 21:57:39 UTC (rev 6355)
@@ -136,9 +136,13 @@
 
     /* Test for singularity */
     if (pivmax < 0.0) {
+#if SCIPY_SPECIFIC_FIX
+        ABORT("[0]: matrix is singular");
+#else
 	fprintf(stderr, "[0]: jcol=%d, SINGULAR!!!\n", jcol);
 	fflush(stderr);
 	exit(1);
+#endif
     }
     if ( pivmax == 0.0 ) {
 	if (diag != EMPTY)
@@ -151,9 +155,13 @@
 	    for (icol = jcol; icol < n; icol++)
 		if (marker[swap[icol]] <= jcol) break;
 	    if (icol >= n) {
+#if SCIPY_SPECIFIC_FIX
+                ABORT("[1]: matrix is singular");
+#else
 		fprintf(stderr, "[1]: jcol=%d, SINGULAR!!!\n", jcol);
 		fflush(stderr);
 		exit(1);
+#endif
 	    }
 
 	    *pivrow = swap[icol];

Modified: trunk/scipy/sparse/linalg/dsolve/SuperLU/SRC/ilu_dpivotL.c
===================================================================
--- trunk/scipy/sparse/linalg/dsolve/SuperLU/SRC/ilu_dpivotL.c	2010-04-27 21:57:15 UTC (rev 6354)
+++ trunk/scipy/sparse/linalg/dsolve/SuperLU/SRC/ilu_dpivotL.c	2010-04-27 21:57:39 UTC (rev 6355)
@@ -134,9 +134,13 @@
 
     /* Test for singularity */
     if (pivmax < 0.0) {
+#if SCIPY_SPECIFIC_FIX
+        ABORT("[0]: matrix is singular");
+#else
 	fprintf(stderr, "[0]: jcol=%d, SINGULAR!!!\n", jcol);
 	fflush(stderr);
 	exit(1);
+#endif
     }
     if ( pivmax == 0.0 ) {
 	if (diag != EMPTY)
@@ -149,9 +153,13 @@
 	    for (icol = jcol; icol < n; icol++)
 		if (marker[swap[icol]] <= jcol) break;
 	    if (icol >= n) {
+#if SCIPY_SPECIFIC_FIX
+                ABORT("[1]: matrix is singular");
+#else
 		fprintf(stderr, "[1]: jcol=%d, SINGULAR!!!\n", jcol);
 		fflush(stderr);
 		exit(1);
+#endif
 	    }
 
 	    *pivrow = swap[icol];

Modified: trunk/scipy/sparse/linalg/dsolve/SuperLU/SRC/ilu_spivotL.c
===================================================================
--- trunk/scipy/sparse/linalg/dsolve/SuperLU/SRC/ilu_spivotL.c	2010-04-27 21:57:15 UTC (rev 6354)
+++ trunk/scipy/sparse/linalg/dsolve/SuperLU/SRC/ilu_spivotL.c	2010-04-27 21:57:39 UTC (rev 6355)
@@ -134,9 +134,13 @@
 
     /* Test for singularity */
     if (pivmax < 0.0) {
+#if SCIPY_SPECIFIC_FIX
+        ABORT("[0]: matrix is singular");
+#else
 	fprintf(stderr, "[0]: jcol=%d, SINGULAR!!!\n", jcol);
 	fflush(stderr);
 	exit(1);
+#endif
     }
     if ( pivmax == 0.0 ) {
 	if (diag != EMPTY)
@@ -149,9 +153,13 @@
 	    for (icol = jcol; icol < n; icol++)
 		if (marker[swap[icol]] <= jcol) break;
 	    if (icol >= n) {
+#if SCIPY_SPECIFIC_FIX
+                ABORT("[1]: matrix is singular");
+#else
 		fprintf(stderr, "[1]: jcol=%d, SINGULAR!!!\n", jcol);
 		fflush(stderr);
 		exit(1);
+#endif
 	    }
 
 	    *pivrow = swap[icol];

Modified: trunk/scipy/sparse/linalg/dsolve/SuperLU/SRC/ilu_zpivotL.c
===================================================================
--- trunk/scipy/sparse/linalg/dsolve/SuperLU/SRC/ilu_zpivotL.c	2010-04-27 21:57:15 UTC (rev 6354)
+++ trunk/scipy/sparse/linalg/dsolve/SuperLU/SRC/ilu_zpivotL.c	2010-04-27 21:57:39 UTC (rev 6355)
@@ -136,9 +136,13 @@
 
     /* Test for singularity */
     if (pivmax < 0.0) {
+#if SCIPY_SPECIFIC_FIX
+        ABORT("[0]: matrix is singular");
+#else
 	fprintf(stderr, "[0]: jcol=%d, SINGULAR!!!\n", jcol);
 	fflush(stderr);
 	exit(1);
+#endif
     }
     if ( pivmax == 0.0 ) {
 	if (diag != EMPTY)
@@ -151,9 +155,13 @@
 	    for (icol = jcol; icol < n; icol++)
 		if (marker[swap[icol]] <= jcol) break;
 	    if (icol >= n) {
+#if SCIPY_SPECIFIC_FIX
+                ABORT("[1]: matrix is singular");
+#else
 		fprintf(stderr, "[1]: jcol=%d, SINGULAR!!!\n", jcol);
 		fflush(stderr);
 		exit(1);
+#endif
 	    }
 
 	    *pivrow = swap[icol];

Modified: trunk/scipy/sparse/linalg/dsolve/SuperLU/scipychanges.txt
===================================================================
--- trunk/scipy/sparse/linalg/dsolve/SuperLU/scipychanges.txt	2010-04-27 21:57:15 UTC (rev 6354)
+++ trunk/scipy/sparse/linalg/dsolve/SuperLU/scipychanges.txt	2010-04-27 21:57:39 UTC (rev 6355)
@@ -17,3 +17,5 @@
 3) BUG: scipy.sparse.dsolve/SuperLU: sprinkle volatile into dlamc/slamc implementation to avoid an infinite loop
 
 4) BUG: sparse.linalg.dsolve/SuperLU: patch SuperLU sources to eliminate a crash for singular matrices for which pivoting fails
+
+5) BUG: parse.linalg.dsolve/SuperLU: patch SuperLU sources to not exit(1) when ILU decomposition encounters singularity; instead, raise a Python exception

Modified: trunk/scipy/sparse/linalg/dsolve/_superluobject.c
===================================================================
--- trunk/scipy/sparse/linalg/dsolve/_superluobject.c	2010-04-27 21:57:15 UTC (rev 6354)
+++ trunk/scipy/sparse/linalg/dsolve/_superluobject.c	2010-04-27 21:57:39 UTC (rev 6355)
@@ -120,8 +120,12 @@
 {
   SUPERLU_FREE(self->perm_r);
   SUPERLU_FREE(self->perm_c);
-  Destroy_SuperNode_Matrix(&self->L);
-  Destroy_CompCol_Matrix(&self->U);
+  if (self->L.Store != NULL) {
+      Destroy_SuperNode_Matrix(&self->L);
+  }
+  if (self->U.Store != NULL) {
+      Destroy_CompCol_Matrix(&self->U);
+  }
   PyObject_Del(self);
 }
 
@@ -291,6 +295,7 @@
   superlu_options_t options;
   SuperLUStat_t stat;
   int panel_size, relax;
+  int trf_finished = 0;
 
   n = A->ncol;
 
@@ -338,6 +343,7 @@
             etree, NULL, lwork, self->perm_c, self->perm_r,
             &self->L, &self->U, &stat, &info);
   }
+  trf_finished = 1;
 
   if (info) {
     if (info < 0)
@@ -360,6 +366,12 @@
   return (PyObject *)self;
 
 fail:
+  if (!trf_finished) {
+      /* Avoid trying to free partially initialized matrices;
+         might leak some memory, but avoids a crash */
+      self->L.Store = NULL;
+      self->U.Store = NULL;
+  }
   SUPERLU_FREE(etree);
   Destroy_CompCol_Permuted(&AC);
   StatFree(&stat);




More information about the Scipy-svn mailing list