Mercurial > hg > octave-lyh

diff liboctave/SparseCmplxQR.cc @ 5681:233d98d95659

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
[project @ 2006-03-16 17:48:55 by dbateman]
author dbateman
date Thu, 16 Mar 2006 17:48:56 +0000
parents 69a4f320d95a
children ace8d8d26933
line wrap: on
line diff
--- a/liboctave/SparseCmplxQR.cc
+++ b/liboctave/SparseCmplxQR.cc
@@ -34,6 +34,8 @@
   OCTAVE_LOCAL_BUFFER (double, buf ## tmp, (2 * (n))); \
   double _Complex *buf = reinterpret_cast<double _Complex *> (buf ## tmp);
 
+#define OCTAVE_C99_ZERO (0. + 0.iF);
+
 SparseComplexQR::SparseComplexQR_rep::SparseComplexQR_rep 
 (const SparseComplexMatrix& a, int order)
 {
@@ -233,6 +235,8 @@
 	  OCTAVE_QUIT;
 	  for (octave_idx_type j = 0; j < b_nr; j++)
 	    Xx[j] = b.xelem(j,i);
+	  for (octave_idx_type j = nr; j < q.S()->m2; j++)
+	    buf[j] = OCTAVE_C99_ZERO;
 	  BEGIN_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
 	  CXSPARSE_ZNAME (_ipvec) 
 	    (nr, q.S()->Pinv, reinterpret_cast<double _Complex *>(Xx), buf);
@@ -262,7 +266,8 @@
       x.resize(nc, b_nc);
       double _Complex *vec = reinterpret_cast<double _Complex *>
 	(x.fortran_vec());
-      OCTAVE_C99_COMPLEX (buf, nc > q.S()->m2 ? nc : q.S()->m2);
+      volatile octave_idx_type nbuf = (nc > q.S()->m2 ? nc : q.S()->m2);
+      OCTAVE_C99_COMPLEX (buf, nbuf);
       OCTAVE_LOCAL_BUFFER (Complex, Xx, b_nr);
       OCTAVE_LOCAL_BUFFER (Complex, B, nr);
       for (octave_idx_type i = 0; i < nr; i++)
@@ -272,6 +277,8 @@
 	  OCTAVE_QUIT;
 	  for (octave_idx_type j = 0; j < b_nr; j++)
 	    Xx[j] = b.xelem(j,i);
+	  for (octave_idx_type j = nr; j < nbuf; j++)
+	    buf[j] = OCTAVE_C99_ZERO;
 	  BEGIN_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
 	  CXSPARSE_ZNAME (_pvec)
 	    (nr, q.S()->Q, reinterpret_cast<double _Complex *>(Xx), buf);
@@ -332,6 +339,8 @@
 	  OCTAVE_QUIT;
 	  for (octave_idx_type j = 0; j < b_nr; j++)
 	    Xx[j] = b.xelem(j,i);
+	  for (octave_idx_type j = nr; j < q.S()->m2; j++)
+	    buf[j] = OCTAVE_C99_ZERO;
 	  BEGIN_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
 	  CXSPARSE_ZNAME (_ipvec) 
 	    (nr, q.S()->Pinv, reinterpret_cast<double _Complex *>(Xx), buf);
@@ -382,8 +391,9 @@
       x.xcidx(0) = 0;
       x_nz = b.nzmax();
       ii = 0;
+      volatile octave_idx_type nbuf = (nc > q.S()->m2 ? nc : q.S()->m2);
       OCTAVE_LOCAL_BUFFER (Complex, Xx, (b_nr > nc ? b_nr : nc));
-      OCTAVE_C99_COMPLEX (buf, nc > q.S()->m2 ? nc : q.S()->m2);
+      OCTAVE_C99_COMPLEX (buf, nbuf);
       OCTAVE_LOCAL_BUFFER (Complex, B, nr);
       for (octave_idx_type i = 0; i < nr; i++)
 	B[i] = conj (reinterpret_cast<Complex *>(q.N()->B) [i]);
@@ -392,6 +402,8 @@
 	  OCTAVE_QUIT;
 	  for (octave_idx_type j = 0; j < b_nr; j++)
 	    Xx[j] = b.xelem(j,i);
+	  for (octave_idx_type j = nr; j < nbuf; j++)
+	    buf[j] = OCTAVE_C99_ZERO;
 	  BEGIN_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
 	  CXSPARSE_ZNAME (_pvec)
 	    (nr, q.S()->Q, reinterpret_cast<double _Complex *>(Xx), buf);
@@ -470,6 +482,8 @@
 	   i++, idx+=nc, bidx+=b_nr)
 	{
 	  OCTAVE_QUIT;
+	  for (octave_idx_type j = nr; j < q.S()->m2; j++)
+	    buf[j] = OCTAVE_C99_ZERO;
 	  BEGIN_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
 	  CXSPARSE_ZNAME (_ipvec) (nr, q.S()->Pinv, bvec + bidx, buf);
 	  END_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
@@ -498,7 +512,8 @@
       x.resize(nc, b_nc);
       double _Complex *vec = reinterpret_cast<double _Complex *>
 	(x.fortran_vec());
-      OCTAVE_C99_COMPLEX (buf, nc > q.S()->m2 ? nc : q.S()->m2);
+      volatile octave_idx_type nbuf = (nc > q.S()->m2 ? nc : q.S()->m2);
+      OCTAVE_C99_COMPLEX (buf, nbuf);
       OCTAVE_LOCAL_BUFFER (Complex, B, nr);
       for (octave_idx_type i = 0; i < nr; i++)
 	B[i] = conj (reinterpret_cast<Complex *>(q.N()->B) [i]);
@@ -506,6 +521,8 @@
 	   i++, idx+=nc, bidx+=b_nr)
 	{
 	  OCTAVE_QUIT;
+	  for (octave_idx_type j = nr; j < nbuf; j++)
+	    buf[j] = OCTAVE_C99_ZERO;
 	  BEGIN_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
 	  CXSPARSE_ZNAME (_pvec) (nr, q.S()->Q, bvec + bidx, buf);
 	  CXSPARSE_ZNAME (_utsolve) (q.N()->U, buf);
@@ -564,6 +581,8 @@
 	  OCTAVE_QUIT;
 	  for (octave_idx_type j = 0; j < b_nr; j++)
 	    Xx[j] = b.xelem(j,i);
+	  for (octave_idx_type j = nr; j < q.S()->m2; j++)
+	    buf[j] = OCTAVE_C99_ZERO;
 	  BEGIN_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
 	  CXSPARSE_ZNAME (_ipvec) 
 	    (nr, q.S()->Pinv, reinterpret_cast<double _Complex *>(Xx), buf);
@@ -614,8 +633,9 @@
       x.xcidx(0) = 0;
       x_nz = b.nzmax();
       ii = 0;
+      volatile octave_idx_type nbuf = (nc > q.S()->m2 ? nc : q.S()->m2);
       OCTAVE_LOCAL_BUFFER (Complex, Xx, (b_nr > nc ? b_nr : nc));
-      OCTAVE_C99_COMPLEX (buf, nc > q.S()->m2 ? nc : q.S()->m2);
+      OCTAVE_C99_COMPLEX (buf, nbuf);
       OCTAVE_LOCAL_BUFFER (Complex, B, nr);
       for (octave_idx_type i = 0; i < nr; i++)
 	B[i] = conj (reinterpret_cast<Complex *>(q.N()->B) [i]);
@@ -624,6 +644,8 @@
 	  OCTAVE_QUIT;
 	  for (octave_idx_type j = 0; j < b_nr; j++)
 	    Xx[j] = b.xelem(j,i);
+	  for (octave_idx_type j = nr; j < nbuf; j++)
+	    buf[j] = OCTAVE_C99_ZERO;
 	  BEGIN_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
 	  CXSPARSE_ZNAME (_pvec)
 	    (nr, q.S()->Q, reinterpret_cast<double _Complex *>(Xx), buf);
@@ -675,4 +697,3 @@
 ;;; mode: C++ ***
 ;;; End: ***
 */
-