commons-commits mailing list archives

Site index · List index
Message view « Date » · « Thread »
Top « Date » · « Thread »
From l..@apache.org
Subject svn commit: r699850 - in /commons/proper/math/branches/MATH_2_0/src: java/org/apache/commons/math/linear/QRDecomposition.java java/org/apache/commons/math/linear/QRDecompositionImpl.java test/org/apache/commons/math/linear/QRDecompositionImplTest.java
Date Sun, 28 Sep 2008 16:17:42 GMT
Author: luc
Date: Sun Sep 28 09:17:41 2008
New Revision: 699850

URL: http://svn.apache.org/viewvc?rev=699850&view=rev
Log:
added a getQT() method to QRDecomposition interface
compute Q by transposing QT rather than the other way round for efficiency
JIRA: MATH-223

Modified:
    commons/proper/math/branches/MATH_2_0/src/java/org/apache/commons/math/linear/QRDecomposition.java
    commons/proper/math/branches/MATH_2_0/src/java/org/apache/commons/math/linear/QRDecompositionImpl.java
    commons/proper/math/branches/MATH_2_0/src/test/org/apache/commons/math/linear/QRDecompositionImplTest.java

Modified: commons/proper/math/branches/MATH_2_0/src/java/org/apache/commons/math/linear/QRDecomposition.java
URL: http://svn.apache.org/viewvc/commons/proper/math/branches/MATH_2_0/src/java/org/apache/commons/math/linear/QRDecomposition.java?rev=699850&r1=699849&r2=699850&view=diff
==============================================================================
--- commons/proper/math/branches/MATH_2_0/src/java/org/apache/commons/math/linear/QRDecomposition.java
(original)
+++ commons/proper/math/branches/MATH_2_0/src/java/org/apache/commons/math/linear/QRDecomposition.java
Sun Sep 28 09:17:41 2008
@@ -58,6 +58,15 @@
     RealMatrix getQ() throws IllegalStateException;
 
     /**
+     * Returns the transpose of the matrix Q of the decomposition.
+     * <p>Q is an orthogonal matrix</p>
+     * @return the Q matrix
+     * @exception IllegalStateException if {@link
+     * DecompositionSolver#decompose(RealMatrix) decompose} has not been called
+     */
+    RealMatrix getQT() throws IllegalStateException;
+
+    /**
      * Returns the Householder reflector vectors.
      * <p>H is a lower trapezoidal matrix whose columns represent
      * each successive Householder reflector vector. This matrix is used

Modified: commons/proper/math/branches/MATH_2_0/src/java/org/apache/commons/math/linear/QRDecompositionImpl.java
URL: http://svn.apache.org/viewvc/commons/proper/math/branches/MATH_2_0/src/java/org/apache/commons/math/linear/QRDecompositionImpl.java?rev=699850&r1=699849&r2=699850&view=diff
==============================================================================
--- commons/proper/math/branches/MATH_2_0/src/java/org/apache/commons/math/linear/QRDecompositionImpl.java
(original)
+++ commons/proper/math/branches/MATH_2_0/src/java/org/apache/commons/math/linear/QRDecompositionImpl.java
Sun Sep 28 09:17:41 2008
@@ -52,6 +52,9 @@
     /** Cached value of Q. */
     private RealMatrix cachedQ;
 
+    /** Cached value of QT. */
+    private RealMatrix cachedQT;
+
     /** Cached value of R. */
     private RealMatrix cachedR;
 
@@ -87,9 +90,10 @@
         final int n = matrix.getColumnDimension();
         qrt = matrix.transpose().getData();
         rDiag = new double[n];
-        cachedQ = null;
-        cachedR = null;
-        cachedH = null;
+        cachedQ  = null;
+        cachedQT = null;
+        cachedR  = null;
+        cachedH  = null;
 
         /*
          * The QR decomposition of a matrix A is calculated using Householder
@@ -191,15 +195,24 @@
     /** {@inheritDoc} */
     public RealMatrix getQ()
         throws IllegalStateException {
-
         if (cachedQ == null) {
+            cachedQ = getQT().transpose();
+        }
+        return cachedQ;
+    }
+
+    /** {@inheritDoc} */
+    public RealMatrix getQT()
+        throws IllegalStateException {
+
+            if (cachedQ == null) {
 
             checkDecomposed();
 
-            // Q is supposed to be m x m
+            // QT is supposed to be m x m
             final int n = qrt.length;
             final int m = qrt[0].length;
-            double[][] q = new double[m][m];
+            double[][] qT = new double[m][m];
 
             /* 
              * Q = Q1 Q2 ... Q_m, so Q is formed by first constructing Q_m and then 
@@ -207,34 +220,35 @@
              * succession to the result 
              */ 
             for (int minor = m - 1; minor >= Math.min(m, n); minor--) {
-                q[minor][minor]=1;
+                qT[minor][minor]=1;
             }
 
             for (int minor = Math.min(m,n)-1; minor >= 0; minor--){
                 final double[] qrtMinor = qrt[minor];
-                q[minor][minor] = 1;
+                qT[minor][minor] = 1;
                 if (qrtMinor[minor] != 0.0) {
                     for (int col = minor; col < m; col++) {
+                        final double[] qTCol = qT[col];
                         double alpha = 0;
                         for (int row = minor; row < m; row++) {
-                            alpha -= q[row][col] * qrtMinor[row];
+                            alpha -= qTCol[row] * qrtMinor[row];
                         }
                         alpha /= rDiag[minor] * qrtMinor[minor];
 
                         for (int row = minor; row < m; row++) {
-                            q[row][col] -= alpha * qrtMinor[row];
+                            qTCol[row] -= alpha * qrtMinor[row];
                         }
                     }
                 }
             }
 
             // cache the matrix for subsequent calls
-            cachedQ = new RealMatrixImpl(q, false);
+            cachedQT = new RealMatrixImpl(qT, false);
 
         }
 
         // return the cached matrix
-        return cachedQ;
+        return cachedQT;
 
     }
 

Modified: commons/proper/math/branches/MATH_2_0/src/test/org/apache/commons/math/linear/QRDecompositionImplTest.java
URL: http://svn.apache.org/viewvc/commons/proper/math/branches/MATH_2_0/src/test/org/apache/commons/math/linear/QRDecompositionImplTest.java?rev=699850&r1=699849&r2=699850&view=diff
==============================================================================
--- commons/proper/math/branches/MATH_2_0/src/test/org/apache/commons/math/linear/QRDecompositionImplTest.java
(original)
+++ commons/proper/math/branches/MATH_2_0/src/test/org/apache/commons/math/linear/QRDecompositionImplTest.java
Sun Sep 28 09:17:41 2008
@@ -109,31 +109,31 @@
     /** test the orthogonality of Q */
     public void testQOrthogonal() {
         RealMatrix matrix = new RealMatrixImpl(testData3x3NonSingular, false);
-        matrix = new QRDecompositionImpl(matrix).getQ();
+        RealMatrix q  = new QRDecompositionImpl(matrix).getQ();
+        RealMatrix qT = new QRDecompositionImpl(matrix).getQT();
         RealMatrix eye = MatrixUtils.createRealIdentityMatrix(3);
-        double norm = matrix.transpose().multiply(matrix).subtract(eye)
-                .getNorm();
+        double norm = qT.multiply(q).subtract(eye).getNorm();
         assertEquals("3x3 nonsingular Q'Q = I", 0, norm, normTolerance);
 
         matrix = new RealMatrixImpl(testData3x3Singular, false);
-        matrix = new QRDecompositionImpl(matrix).getQ();
+        q  = new QRDecompositionImpl(matrix).getQ();
+        qT = new QRDecompositionImpl(matrix).getQT();
         eye = MatrixUtils.createRealIdentityMatrix(3);
-        norm = matrix.transpose().multiply(matrix).subtract(eye)
-                .getNorm();
+        norm = qT.multiply(q).subtract(eye).getNorm();
         assertEquals("3x3 singular Q'Q = I", 0, norm, normTolerance);
 
         matrix = new RealMatrixImpl(testData3x4, false);
-        matrix = new QRDecompositionImpl(matrix).getQ();
+        q  = new QRDecompositionImpl(matrix).getQ();
+        qT = new QRDecompositionImpl(matrix).getQT();
         eye = MatrixUtils.createRealIdentityMatrix(3);
-        norm = matrix.transpose().multiply(matrix).subtract(eye)
-                .getNorm();
+        norm = qT.multiply(q).subtract(eye).getNorm();
         assertEquals("3x4 Q'Q = I", 0, norm, normTolerance);
 
         matrix = new RealMatrixImpl(testData4x3, false);
-        matrix = new QRDecompositionImpl(matrix).getQ();
+        q  = new QRDecompositionImpl(matrix).getQ();
+        qT = new QRDecompositionImpl(matrix).getQT();
         eye = MatrixUtils.createRealIdentityMatrix(4);
-        norm = matrix.transpose().multiply(matrix).subtract(eye)
-                .getNorm();
+        norm = qT.multiply(q).subtract(eye).getNorm();
         assertEquals("4x3 Q'Q = I", 0, norm, normTolerance);
     }
 
@@ -337,6 +337,8 @@
         // check values against known references
         RealMatrix q = qr.getQ();
         assertEquals(0, q.subtract(qRef).getNorm(), 1.0e-13);
+        RealMatrix qT = qr.getQT();
+        assertEquals(0, qT.subtract(qRef.transpose()).getNorm(), 1.0e-13);
         RealMatrix r = qr.getR();
         assertEquals(0, r.subtract(rRef).getNorm(), 1.0e-13);
         RealMatrix h = qr.getH();



Mime
View raw message