diff liboctave/dbleQR.cc @ 539:5ec10a984241

[project @ 1994-07-21 22:40:04 by jwe]
author jwe
date Thu, 21 Jul 1994 22:42:09 +0000
parents 3d4b4f0fa5ba
children 714fd17fca28
line wrap: on
line diff
--- a/liboctave/dbleQR.cc	Thu Jul 21 22:40:04 1994 +0000
+++ b/liboctave/dbleQR.cc	Thu Jul 21 22:42:09 1994 +0000
@@ -41,9 +41,12 @@
 
   int F77_FCN (dorgqr) (const int*, const int*, const int*, double*,
 			const int*, double*, double*, const int*, int*);
+
+  int F77_FCN (dgeqpf) (const int*, const int*, double*, const int*,
+			int*, double*, double*, int*);
 }
 
-QR::QR (const Matrix& a)
+QR::QR (const Matrix& a, QR::type qr_type)
 {
   int m = a.rows ();
   int n = a.cols ();
@@ -73,23 +76,47 @@
 
   delete [] work;
 
-  r.resize (m, n, 0.0);
-  for (int j = 0; j < n; j++)
+  if (qr_type == QR::raw)
+    {
+      for (int j = 0; j < min_mn; j++)
+	{
+	  int limit = j < min_mn - 1 ? j : min_mn - 1;
+	  for (int i = limit + 1; i < m; i++)
+	    tmp_data[m*j+i] *= tau[j];
+	}
+    }
+  else
     {
-      int limit = j < min_mn-1 ? j : min_mn-1;
-      for (int i = 0; i <= limit; i++)
-	r.elem (i, j) = tmp_data[m*j+i];
-    }
+      int m2;
+      if (qr_type == QR::economy && m > n)
+	{
+	  m2 = n;
+	  r.resize (n, n, 0.0);
+	}
+      else
+	{
+	  m2 = m;
+	  r.resize (m, n, 0.0);
+	}
 
-  lwork = 32*m;
-  work = new double[lwork];
-
-  F77_FCN (dorgqr) (&m, &m, &min_mn, tmp_data, &m, tau, work, &lwork, &info);
+      for (int j = 0; j < n; j++)
+	{
+	  int limit = j < min_mn-1 ? j : min_mn-1;
+	  for (int i = 0; i <= limit; i++)
+	    r.elem (i, j) = tmp_data[m*j+i];
+	}
 
-  q = Matrix (tmp_data, m, m);
+      lwork = 32*m;
+      work = new double[lwork];
+
+      F77_FCN (dorgqr) (&m, &m, &min_mn, tmp_data, &m, tau, work,
+			&lwork, &info);
 
-  delete [] tau;
-  delete [] work;
+      q = Matrix (tmp_data, m, m);
+
+      delete [] tau;
+      delete [] work;
+    }
 }
 
 /*