View Javadoc
1   /*
2    * Licensed to the Apache Software Foundation (ASF) under one or more
3    * contributor license agreements.  See the NOTICE file distributed with
4    * this work for additional information regarding copyright ownership.
5    * The ASF licenses this file to You under the Apache License, Version 2.0
6    * (the "License"); you may not use this file except in compliance with
7    * the License.  You may obtain a copy of the License at
8    *
9    *      https://www.apache.org/licenses/LICENSE-2.0
10   *
11   * Unless required by applicable law or agreed to in writing, software
12   * distributed under the License is distributed on an "AS IS" BASIS,
13   * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
14   * See the License for the specific language governing permissions and
15   * limitations under the License.
16   */
17  
18  /*
19   * This is not the original file distributed by the Apache Software Foundation
20   * It has been modified by the Hipparchus project
21   */
22  
23  package org.hipparchus.linear;
24  
25  import java.util.Random;
26  
27  import org.hipparchus.exception.MathIllegalArgumentException;
28  import org.junit.Assert;
29  import org.junit.Test;
30  
31  
32  public class RRQRDecompositionTest {
33      private double[][] testData3x3NonSingular = {
34              { 12, -51, 4 },
35              { 6, 167, -68 },
36              { -4, 24, -41 }, };
37  
38      private double[][] testData3x3Singular = {
39              { 1, 4, 7, },
40              { 2, 5, 8, },
41              { 3, 6, 9, }, };
42  
43      private double[][] testData3x4 = {
44              { 12, -51, 4, 1 },
45              { 6, 167, -68, 2 },
46              { -4, 24, -41, 3 }, };
47  
48      private double[][] testData4x3 = {
49              { 12, -51, 4, },
50              { 6, 167, -68, },
51              { -4, 24, -41, },
52              { -5, 34, 7, }, };
53  
54      private static final double entryTolerance = 10e-16;
55  
56      private static final double normTolerance = 10e-14;
57  
58      /** test dimensions */
59      @Test
60      public void testDimensions() {
61          checkDimension(MatrixUtils.createRealMatrix(testData3x3NonSingular));
62  
63          checkDimension(MatrixUtils.createRealMatrix(testData4x3));
64  
65          checkDimension(MatrixUtils.createRealMatrix(testData3x4));
66  
67          Random r = new Random(643895747384642l);
68          int    p = (5 * BlockRealMatrix.BLOCK_SIZE) / 4;
69          int    q = (7 * BlockRealMatrix.BLOCK_SIZE) / 4;
70          checkDimension(createTestMatrix(r, p, q));
71          checkDimension(createTestMatrix(r, q, p));
72  
73      }
74  
75      private void checkDimension(RealMatrix m) {
76          int rows = m.getRowDimension();
77          int columns = m.getColumnDimension();
78          RRQRDecomposition qr = new RRQRDecomposition(m);
79          Assert.assertEquals(rows,    qr.getQ().getRowDimension());
80          Assert.assertEquals(rows,    qr.getQ().getColumnDimension());
81          Assert.assertEquals(rows,    qr.getR().getRowDimension());
82          Assert.assertEquals(columns, qr.getR().getColumnDimension());
83      }
84  
85      /** test AP = QR */
86      @Test
87      public void testAPEqualQR() {
88          checkAPEqualQR(MatrixUtils.createRealMatrix(testData3x3NonSingular));
89  
90          checkAPEqualQR(MatrixUtils.createRealMatrix(testData3x3Singular));
91  
92          checkAPEqualQR(MatrixUtils.createRealMatrix(testData3x4));
93  
94          checkAPEqualQR(MatrixUtils.createRealMatrix(testData4x3));
95  
96          Random r = new Random(643895747384642l);
97          int    p = (5 * BlockRealMatrix.BLOCK_SIZE) / 4;
98          int    q = (7 * BlockRealMatrix.BLOCK_SIZE) / 4;
99          checkAPEqualQR(createTestMatrix(r, p, q));
100 
101         checkAPEqualQR(createTestMatrix(r, q, p));
102 
103     }
104 
105     private void checkAPEqualQR(RealMatrix m) {
106         RRQRDecomposition rrqr = new RRQRDecomposition(m);
107         double norm = rrqr.getQ().multiply(rrqr.getR()).subtract(m.multiply(rrqr.getP())).getNorm1();
108         Assert.assertEquals(0, norm, normTolerance);
109     }
110 
111     /** test the orthogonality of Q */
112     @Test
113     public void testQOrthogonal() {
114         checkQOrthogonal(MatrixUtils.createRealMatrix(testData3x3NonSingular));
115 
116         checkQOrthogonal(MatrixUtils.createRealMatrix(testData3x3Singular));
117 
118         checkQOrthogonal(MatrixUtils.createRealMatrix(testData3x4));
119 
120         checkQOrthogonal(MatrixUtils.createRealMatrix(testData4x3));
121 
122         Random r = new Random(643895747384642l);
123         int    p = (5 * BlockRealMatrix.BLOCK_SIZE) / 4;
124         int    q = (7 * BlockRealMatrix.BLOCK_SIZE) / 4;
125         checkQOrthogonal(createTestMatrix(r, p, q));
126 
127         checkQOrthogonal(createTestMatrix(r, q, p));
128 
129     }
130 
131     private void checkQOrthogonal(RealMatrix m) {
132         RRQRDecomposition qr = new RRQRDecomposition(m);
133         RealMatrix eye = MatrixUtils.createRealIdentityMatrix(m.getRowDimension());
134         double norm = qr.getQT().multiply(qr.getQ()).subtract(eye).getNorm1();
135         Assert.assertEquals(0, norm, normTolerance);
136     }
137 
138     /** test that R is upper triangular */
139     @Test
140     public void testRUpperTriangular() {
141         RealMatrix matrix = MatrixUtils.createRealMatrix(testData3x3NonSingular);
142         checkUpperTriangular(new RRQRDecomposition(matrix).getR());
143 
144         matrix = MatrixUtils.createRealMatrix(testData3x3Singular);
145         checkUpperTriangular(new RRQRDecomposition(matrix).getR());
146 
147         matrix = MatrixUtils.createRealMatrix(testData3x4);
148         checkUpperTriangular(new RRQRDecomposition(matrix).getR());
149 
150         matrix = MatrixUtils.createRealMatrix(testData4x3);
151         checkUpperTriangular(new RRQRDecomposition(matrix).getR());
152 
153         Random r = new Random(643895747384642l);
154         int    p = (5 * BlockRealMatrix.BLOCK_SIZE) / 4;
155         int    q = (7 * BlockRealMatrix.BLOCK_SIZE) / 4;
156         matrix = createTestMatrix(r, p, q);
157         checkUpperTriangular(new RRQRDecomposition(matrix).getR());
158 
159         matrix = createTestMatrix(r, p, q);
160         checkUpperTriangular(new RRQRDecomposition(matrix).getR());
161 
162     }
163 
164     private void checkUpperTriangular(RealMatrix m) {
165         m.walkInOptimizedOrder(new DefaultRealMatrixPreservingVisitor() {
166             @Override
167             public void visit(int row, int column, double value) {
168                 if (column < row) {
169                     Assert.assertEquals(0.0, value, entryTolerance);
170                 }
171             }
172         });
173     }
174 
175     /** test that H is trapezoidal */
176     @Test
177     public void testHTrapezoidal() {
178         RealMatrix matrix = MatrixUtils.createRealMatrix(testData3x3NonSingular);
179         checkTrapezoidal(new RRQRDecomposition(matrix).getH());
180 
181         matrix = MatrixUtils.createRealMatrix(testData3x3Singular);
182         checkTrapezoidal(new RRQRDecomposition(matrix).getH());
183 
184         matrix = MatrixUtils.createRealMatrix(testData3x4);
185         checkTrapezoidal(new RRQRDecomposition(matrix).getH());
186 
187         matrix = MatrixUtils.createRealMatrix(testData4x3);
188         checkTrapezoidal(new RRQRDecomposition(matrix).getH());
189 
190         Random r = new Random(643895747384642l);
191         int    p = (5 * BlockRealMatrix.BLOCK_SIZE) / 4;
192         int    q = (7 * BlockRealMatrix.BLOCK_SIZE) / 4;
193         matrix = createTestMatrix(r, p, q);
194         checkTrapezoidal(new RRQRDecomposition(matrix).getH());
195 
196         matrix = createTestMatrix(r, p, q);
197         checkTrapezoidal(new RRQRDecomposition(matrix).getH());
198 
199     }
200 
201     private void checkTrapezoidal(RealMatrix m) {
202         m.walkInOptimizedOrder(new DefaultRealMatrixPreservingVisitor() {
203             @Override
204             public void visit(int row, int column, double value) {
205                 if (column > row) {
206                     Assert.assertEquals(0.0, value, entryTolerance);
207                 }
208             }
209         });
210     }
211 
212     @Test(expected=MathIllegalArgumentException.class)
213     public void testNonInvertible() {
214         RRQRDecomposition qr =
215             new RRQRDecomposition(MatrixUtils.createRealMatrix(testData3x3Singular), 3.0e-16);
216         qr.getSolver().getInverse();
217     }
218 
219     private RealMatrix createTestMatrix(final Random r, final int rows, final int columns) {
220         RealMatrix m = MatrixUtils.createRealMatrix(rows, columns);
221         m.walkInOptimizedOrder(new DefaultRealMatrixChangingVisitor(){
222             @Override
223             public double visit(int row, int column, double value) {
224                 return 2.0 * r.nextDouble() - 1.0;
225             }
226         });
227         return m;
228     }
229 
230     /** test the rank is returned correctly */
231     @Test
232     public void testRank() {
233         double[][] d = { { 1, 1, 1 }, { 0, 0, 0 }, { 1, 2, 3 } };
234         RealMatrix m = new Array2DRowRealMatrix(d);
235         RRQRDecomposition qr = new RRQRDecomposition(m);
236         Assert.assertEquals(2, qr.getRank(1.0e-16));
237     }
238 
239 }