[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]
crash in test_sparse
|
From: |
John W. Eaton |
|
Subject: |
crash in test_sparse |
|
Date: |
Tue, 30 Oct 2007 15:25:52 -0400 |
On 30-Oct-2007, John W. Eaton wrote:
| In the last few days I started noticing this problem while running
| make check:
|
| Fixed test scripts:
|
| test_sparse.m ..........................................error: memory
exhausted or requested size too large for range of Octave's index type --
execution of /home/jwe/src/octave-trunk/test/fntests.m failed
| error: assert (pds \ xs,sparse (pdf \ xs, true),100 * eps) expected
| Compressed Column Sparse (rows = 11, cols = 11, nnz = 11)
|
| (2, 1) -> 0.50000
| (1, 2) -> 1
| (3, 3) -> 0.33333
| (4, 4) -> 0.25000
| (5, 5) -> 0.20000
| (6, 6) -> 0.16667
| (7, 7) -> 0.14286
| (8, 8) -> 0.12500
| (9, 9) -> 0.11111
| (10, 10) -> 0.10000
| (11, 11) -> 1
| but got
| Compressed Column Sparse (rows = 10, cols = 11, nnz = 10)
|
| (2, 1) -> 0.50000
| (1, 2) -> 1
| (3, 3) -> 0.33333
| (4, 4) -> 0.25000
| (5, 5) -> 0.20000
| (6, 6) -> 0.16667
| (7, 7) -> 0.14286
| (8, 8) -> 0.12500
| (9, 9) -> 0.11111
| (10, 10) -> 0.10000
| Dimensions don't match
|
|
| Is anyone else seeing this problem?
|
| I do see it with a fresh build of the current CVS sources but not with
| a fresh build fo 2.9.15 so I'm trying to isolate which change caused
| this, but it would be helpful if someone else could confirm that they
| also see the bug.
I think the following patch fixes the problem. It seems that the
interface to ZGELSD has apparently changed. In LAPACK 3.1.1 (the
version we have in libcruft/lapack) ZGELSD returns computed sizes or
RWORK and IWORK but this does not seem to be the way my copy of
ATLAS/LAPACK on Debian currently works. So it seems we need to
computed these sizes ourselves until the rest of the world catches up
to the new LAPACK.
I guess the lesson is that it is not safe to trust the comments in the
lapack headers because the interface may have changed from previous
versions and we can't assume that ATLAS or the vendor LAPACK libraries
are up to date.
jwe
liboctave/ChangeLog:
2007-10-30 John W. Eaton <address@hidden>
* CMatrix.cc (lssolve): Compute size of rwork and iwork arrays.
* dMatrix.cc (lssolve): Compute size of iwork array.
Index: liboctave/CMatrix.cc
===================================================================
RCS file: /cvs/octave/liboctave/CMatrix.cc,v
retrieving revision 1.140
diff -u -u -r1.140 CMatrix.cc
--- liboctave/CMatrix.cc 29 Oct 2007 18:09:57 -0000 1.140
+++ liboctave/CMatrix.cc 30 Oct 2007 19:25:02 -0000
@@ -2491,13 +2491,38 @@
octave_idx_type lwork = -1;
Array<Complex> work (1);
- Array<double> rwork (1);
- Array<octave_idx_type> iwork (1);
+
+ // FIXME: Can SMLSIZ be other than 25?
+ octave_idx_type smlsiz = 25;
+
+ // We compute the size of rwork and iwork because ZGELSD in
+ // older versions of LAPACK does not return them on a query
+ // call.
+#if defined (HAVE_LOG2)
+ double tmp = log2 (minmn) / static_cast<double> (smlsiz+1) + 1;
+#else
+ double tmp = log (minmn) / static_cast<double> (smlsiz+1) / log (2) + 1;
+#endif
+ octave_idx_type nlvl = static_cast<int> (tmp);
+ if (nlvl < 0)
+ nlvl = 0;
+
+ octave_idx_type lrwork = minmn*(10 + 2*smlsiz + 8*nlvl)
+ + 3*smlsiz*nrhs + (smlsiz+1)*(smlsiz+1);
+ if (lrwork < 1)
+ lrwork = 1;
+ Array<double> rwork (lrwork);
+ double *prwork = rwork.fortran_vec ();
+
+ octave_idx_type liwork = 3 * minmn * nlvl + 11 * minmn;
+ if (liwork < 1)
+ liwork = 1;
+ Array<octave_idx_type> iwork (liwork);
+ octave_idx_type* piwork = iwork.fortran_vec ();
F77_XFCN (zgelsd, ZGELSD, (m, n, nrhs, tmp_data, m, pretval, maxmn,
ps, rcond, rank, work.fortran_vec (),
- lwork, rwork.fortran_vec (),
- iwork.fortran_vec (), info));
+ lwork, prwork, piwork, info));
if (f77_exception_encountered)
(*current_liboctave_error_handler)
@@ -2506,14 +2531,11 @@
{
lwork = static_cast<octave_idx_type> (std::real (work(0)));
work.resize (lwork);
- rwork.resize (static_cast<octave_idx_type> (rwork(0)));
- iwork.resize (iwork(0));
F77_XFCN (zgelsd, ZGELSD, (m, n, nrhs, tmp_data, m, pretval,
maxmn, ps, rcond, rank,
work.fortran_vec (), lwork,
- rwork.fortran_vec (),
- iwork.fortran_vec (), info));
+ prwork, piwork, info));
if (f77_exception_encountered)
(*current_liboctave_error_handler)
@@ -2529,6 +2551,8 @@
rcond = 0.0;
else
rcond = s.elem (minmn - 1) / s.elem (0);
+
+ retval.resize (n, nrhs);
}
}
}
@@ -2637,13 +2661,38 @@
octave_idx_type lwork = -1;
Array<Complex> work (1);
- Array<double> rwork (1);
- Array<octave_idx_type> iwork (1);
+
+ // FIXME: Can SMLSIZ be other than 25?
+ octave_idx_type smlsiz = 25;
+
+ // We compute the size of rwork and iwork because ZGELSD in
+ // older versions of LAPACK does not return them on a query
+ // call.
+#if defined (HAVE_LOG2)
+ double tmp = log2 (minmn) / static_cast<double> (smlsiz+1) + 1;
+#else
+ double tmp = log (minmn) / static_cast<double> (smlsiz+1) / log (2) + 1;
+#endif
+ octave_idx_type nlvl = static_cast<int> (tmp);
+ if (nlvl < 0)
+ nlvl = 0;
+
+ octave_idx_type lrwork = minmn*(10 + 2*smlsiz + 8*nlvl)
+ + 3*smlsiz*nrhs + (smlsiz+1)*(smlsiz+1);
+ if (lrwork < 1)
+ lrwork = 1;
+ Array<double> rwork (lrwork);
+ double *prwork = rwork.fortran_vec ();
+
+ octave_idx_type liwork = 3 * minmn * nlvl + 11 * minmn;
+ if (liwork < 1)
+ liwork = 1;
+ Array<octave_idx_type> iwork (liwork);
+ octave_idx_type* piwork = iwork.fortran_vec ();
F77_XFCN (zgelsd, ZGELSD, (m, n, nrhs, tmp_data, m, pretval, maxmn,
ps, rcond, rank, work.fortran_vec (),
- lwork, rwork.fortran_vec (),
- iwork.fortran_vec (), info));
+ lwork, prwork, piwork, info));
if (f77_exception_encountered)
(*current_liboctave_error_handler)
@@ -2658,8 +2707,7 @@
F77_XFCN (zgelsd, ZGELSD, (m, n, nrhs, tmp_data, m, pretval,
maxmn, ps, rcond, rank,
work.fortran_vec (), lwork,
- rwork.fortran_vec (),
- iwork.fortran_vec (), info));
+ prwork, piwork, info));
if (f77_exception_encountered)
(*current_liboctave_error_handler)
@@ -2675,6 +2723,8 @@
rcond = 0.0;
else
rcond = s.elem (minmn - 1) / s.elem (0);
+
+ retval.resize (n, nrhs);
}
}
}
Index: liboctave/dMatrix.cc
===================================================================
RCS file: /cvs/octave/liboctave/dMatrix.cc,v
retrieving revision 1.143
diff -u -u -r1.143 dMatrix.cc
--- liboctave/dMatrix.cc 29 Oct 2007 18:09:57 -0000 1.143
+++ liboctave/dMatrix.cc 30 Oct 2007 19:25:03 -0000
@@ -2104,7 +2104,22 @@
Array<double> work (1);
// FIXME: Can SMLSIZ be other than 25?
- octave_idx_type liwork = 3 * minmn * 25 + 11 * minmn;
+ octave_idx_type smlsiz = 25;
+
+ // We compute the size of iwork because DGELSD in older versions
+ // of LAPACK does not return it on a query call.
+#if defined (HAVE_LOG2)
+ double tmp = log2 (minmn) / static_cast<double> (smlsiz+1) + 1;
+#else
+ double tmp = log (minmn) / static_cast<double> (smlsiz+1) / log (2) + 1;
+#endif
+ octave_idx_type nlvl = static_cast<int> (tmp);
+ if (nlvl < 0)
+ nlvl = 0;
+
+ octave_idx_type liwork = 3 * minmn * nlvl + 11 * minmn;
+ if (liwork < 1)
+ liwork = 1;
Array<octave_idx_type> iwork (liwork);
octave_idx_type* piwork = iwork.fortran_vec ();
@@ -2137,6 +2152,8 @@
rcond = 0.0;
else
rcond = s.elem (minmn - 1) / s.elem (0);
+
+ retval.resize (n, nrhs);
}
}
}
@@ -2250,7 +2267,22 @@
Array<double> work (1);
// FIXME: Can SMLSIZ be other than 25?
- octave_idx_type liwork = 3 * minmn * 25 + 11 * minmn;
+ octave_idx_type smlsiz = 25;
+
+ // We compute the size of iwork because DGELSD in older versions
+ // of LAPACK does not return it on a query call.
+#if defined (HAVE_LOG2)
+ double tmp = log2 (minmn) / static_cast<double> (smlsiz+1) + 1;
+#else
+ double tmp = log (minmn) / static_cast<double> (smlsiz+1) / log (2) + 1;
+#endif
+ octave_idx_type nlvl = static_cast<int> (tmp);
+ if (nlvl < 0)
+ nlvl = 0;
+
+ octave_idx_type liwork = 3 * minmn * nlvl + 11 * minmn;
+ if (liwork < 1)
+ liwork = 1;
Array<octave_idx_type> iwork (liwork);
octave_idx_type* piwork = iwork.fortran_vec ();
@@ -2284,6 +2316,8 @@
else
rcond = s.elem (minmn - 1) / s.elem (0);
}
+
+ retval.resize (n, nrhs);
}
}