diff liboctave/dbleSVD.cc @ 3336:08ad797989f8

[project @ 1999-11-03 21:41:34 by jwe]
author jwe
date Wed, 03 Nov 1999 21:41:35 +0000
parents f39b97e13cf2
children d14c483b3c12
line wrap: on
line diff
--- a/liboctave/dbleSVD.cc	Wed Nov 03 21:05:01 1999 +0000
+++ b/liboctave/dbleSVD.cc	Wed Nov 03 21:41:35 1999 +0000
@@ -80,7 +80,6 @@
   double *tmp_data = atmp.fortran_vec ();
 
   int min_mn = m < n ? m : n;
-  int max_mn = m > n ? m : n;
 
   char jobu = 'A';
   char jobv = 'A';
@@ -131,23 +130,35 @@
 
   double *vt = right_sm.fortran_vec ();
 
-  int tmp1 = 3*min_mn + max_mn;
-  int tmp2 = 5*min_mn - 4;
-  int lwork = tmp1 > tmp2 ? tmp1 : tmp2;
+  // Ask DGESVD what the dimension of WORK should be.
 
-  Array<double> work (lwork);
-  double *pwork = work.fortran_vec ();
+  int lwork = -1;
+
+  Array<double> work (1);
 
   F77_XFCN (dgesvd, DGESVD, (&jobu, &jobv, m, n, tmp_data, m, s_vec,
-			     u, m, vt, nrow_vt, pwork, lwork, info,
-			     1L, 1L));
+			     u, m, vt, nrow_vt, work.fortran_vec (),
+			     lwork, info, 1L, 1L));
 
   if (f77_exception_encountered)
     (*current_liboctave_error_handler) ("unrecoverable error in dgesvd");
   else
     {
-      if (! (jobv == 'N' || jobv == 'O'))
-	right_sm = right_sm.transpose ();
+      lwork = static_cast<int> (work(0));
+      work.resize (lwork);
+
+      F77_XFCN (dgesvd, DGESVD, (&jobu, &jobv, m, n, tmp_data, m,
+				 s_vec, u, m, vt, nrow_vt,
+				 work.fortran_vec (), lwork, info, 1L,
+				 1L));
+
+      if (f77_exception_encountered)
+	(*current_liboctave_error_handler) ("unrecoverable error in dgesvd");
+      else
+	{
+	  if (! (jobv == 'N' || jobv == 'O'))
+	    right_sm = right_sm.transpose ();
+	}
     }
 
   return info;