1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23 package org.hipparchus.linear;
24
25 import org.hipparchus.exception.MathIllegalArgumentException;
26 import org.junit.jupiter.api.Test;
27
28 import java.util.Random;
29
30 import static org.junit.jupiter.api.Assertions.assertEquals;
31 import static org.junit.jupiter.api.Assertions.assertThrows;
32
33
34 class RRQRDecompositionTest {
35 private double[][] testData3x3NonSingular = {
36 { 12, -51, 4 },
37 { 6, 167, -68 },
38 { -4, 24, -41 }, };
39
40 private double[][] testData3x3Singular = {
41 { 1, 4, 7, },
42 { 2, 5, 8, },
43 { 3, 6, 9, }, };
44
45 private double[][] testData3x4 = {
46 { 12, -51, 4, 1 },
47 { 6, 167, -68, 2 },
48 { -4, 24, -41, 3 }, };
49
50 private double[][] testData4x3 = {
51 { 12, -51, 4, },
52 { 6, 167, -68, },
53 { -4, 24, -41, },
54 { -5, 34, 7, }, };
55
56 private static final double entryTolerance = 10e-16;
57
58 private static final double normTolerance = 10e-14;
59
60
61 @Test
62 void testDimensions() {
63 checkDimension(MatrixUtils.createRealMatrix(testData3x3NonSingular));
64
65 checkDimension(MatrixUtils.createRealMatrix(testData4x3));
66
67 checkDimension(MatrixUtils.createRealMatrix(testData3x4));
68
69 Random r = new Random(643895747384642l);
70 int p = (5 * BlockRealMatrix.BLOCK_SIZE) / 4;
71 int q = (7 * BlockRealMatrix.BLOCK_SIZE) / 4;
72 checkDimension(createTestMatrix(r, p, q));
73 checkDimension(createTestMatrix(r, q, p));
74
75 }
76
77 private void checkDimension(RealMatrix m) {
78 int rows = m.getRowDimension();
79 int columns = m.getColumnDimension();
80 RRQRDecomposition qr = new RRQRDecomposition(m);
81 assertEquals(rows, qr.getQ().getRowDimension());
82 assertEquals(rows, qr.getQ().getColumnDimension());
83 assertEquals(rows, qr.getR().getRowDimension());
84 assertEquals(columns, qr.getR().getColumnDimension());
85 }
86
87
88 @Test
89 void testAPEqualQR() {
90 checkAPEqualQR(MatrixUtils.createRealMatrix(testData3x3NonSingular));
91
92 checkAPEqualQR(MatrixUtils.createRealMatrix(testData3x3Singular));
93
94 checkAPEqualQR(MatrixUtils.createRealMatrix(testData3x4));
95
96 checkAPEqualQR(MatrixUtils.createRealMatrix(testData4x3));
97
98 Random r = new Random(643895747384642l);
99 int p = (5 * BlockRealMatrix.BLOCK_SIZE) / 4;
100 int q = (7 * BlockRealMatrix.BLOCK_SIZE) / 4;
101 checkAPEqualQR(createTestMatrix(r, p, q));
102
103 checkAPEqualQR(createTestMatrix(r, q, p));
104
105 }
106
107 private void checkAPEqualQR(RealMatrix m) {
108 RRQRDecomposition rrqr = new RRQRDecomposition(m);
109 double norm = rrqr.getQ().multiply(rrqr.getR()).subtract(m.multiply(rrqr.getP())).getNorm1();
110 assertEquals(0, norm, normTolerance);
111 }
112
113
114 @Test
115 void testQOrthogonal() {
116 checkQOrthogonal(MatrixUtils.createRealMatrix(testData3x3NonSingular));
117
118 checkQOrthogonal(MatrixUtils.createRealMatrix(testData3x3Singular));
119
120 checkQOrthogonal(MatrixUtils.createRealMatrix(testData3x4));
121
122 checkQOrthogonal(MatrixUtils.createRealMatrix(testData4x3));
123
124 Random r = new Random(643895747384642l);
125 int p = (5 * BlockRealMatrix.BLOCK_SIZE) / 4;
126 int q = (7 * BlockRealMatrix.BLOCK_SIZE) / 4;
127 checkQOrthogonal(createTestMatrix(r, p, q));
128
129 checkQOrthogonal(createTestMatrix(r, q, p));
130
131 }
132
133 private void checkQOrthogonal(RealMatrix m) {
134 RRQRDecomposition qr = new RRQRDecomposition(m);
135 RealMatrix eye = MatrixUtils.createRealIdentityMatrix(m.getRowDimension());
136 double norm = qr.getQT().multiply(qr.getQ()).subtract(eye).getNorm1();
137 assertEquals(0, norm, normTolerance);
138 }
139
140
141 @Test
142 void testRUpperTriangular() {
143 RealMatrix matrix = MatrixUtils.createRealMatrix(testData3x3NonSingular);
144 checkUpperTriangular(new RRQRDecomposition(matrix).getR());
145
146 matrix = MatrixUtils.createRealMatrix(testData3x3Singular);
147 checkUpperTriangular(new RRQRDecomposition(matrix).getR());
148
149 matrix = MatrixUtils.createRealMatrix(testData3x4);
150 checkUpperTriangular(new RRQRDecomposition(matrix).getR());
151
152 matrix = MatrixUtils.createRealMatrix(testData4x3);
153 checkUpperTriangular(new RRQRDecomposition(matrix).getR());
154
155 Random r = new Random(643895747384642l);
156 int p = (5 * BlockRealMatrix.BLOCK_SIZE) / 4;
157 int q = (7 * BlockRealMatrix.BLOCK_SIZE) / 4;
158 matrix = createTestMatrix(r, p, q);
159 checkUpperTriangular(new RRQRDecomposition(matrix).getR());
160
161 matrix = createTestMatrix(r, p, q);
162 checkUpperTriangular(new RRQRDecomposition(matrix).getR());
163
164 }
165
166 private void checkUpperTriangular(RealMatrix m) {
167 m.walkInOptimizedOrder(new DefaultRealMatrixPreservingVisitor() {
168 @Override
169 public void visit(int row, int column, double value) {
170 if (column < row) {
171 assertEquals(0.0, value, entryTolerance);
172 }
173 }
174 });
175 }
176
177
178 @Test
179 void testHTrapezoidal() {
180 RealMatrix matrix = MatrixUtils.createRealMatrix(testData3x3NonSingular);
181 checkTrapezoidal(new RRQRDecomposition(matrix).getH());
182
183 matrix = MatrixUtils.createRealMatrix(testData3x3Singular);
184 checkTrapezoidal(new RRQRDecomposition(matrix).getH());
185
186 matrix = MatrixUtils.createRealMatrix(testData3x4);
187 checkTrapezoidal(new RRQRDecomposition(matrix).getH());
188
189 matrix = MatrixUtils.createRealMatrix(testData4x3);
190 checkTrapezoidal(new RRQRDecomposition(matrix).getH());
191
192 Random r = new Random(643895747384642l);
193 int p = (5 * BlockRealMatrix.BLOCK_SIZE) / 4;
194 int q = (7 * BlockRealMatrix.BLOCK_SIZE) / 4;
195 matrix = createTestMatrix(r, p, q);
196 checkTrapezoidal(new RRQRDecomposition(matrix).getH());
197
198 matrix = createTestMatrix(r, p, q);
199 checkTrapezoidal(new RRQRDecomposition(matrix).getH());
200
201 }
202
203 private void checkTrapezoidal(RealMatrix m) {
204 m.walkInOptimizedOrder(new DefaultRealMatrixPreservingVisitor() {
205 @Override
206 public void visit(int row, int column, double value) {
207 if (column > row) {
208 assertEquals(0.0, value, entryTolerance);
209 }
210 }
211 });
212 }
213
214 @Test
215 void testNonInvertible() {
216 assertThrows(MathIllegalArgumentException.class, () -> {
217 RRQRDecomposition qr =
218 new RRQRDecomposition(MatrixUtils.createRealMatrix(testData3x3Singular), 3.0e-16);
219 qr.getSolver().getInverse();
220 });
221 }
222
223 private RealMatrix createTestMatrix(final Random r, final int rows, final int columns) {
224 RealMatrix m = MatrixUtils.createRealMatrix(rows, columns);
225 m.walkInOptimizedOrder(new DefaultRealMatrixChangingVisitor(){
226 @Override
227 public double visit(int row, int column, double value) {
228 return 2.0 * r.nextDouble() - 1.0;
229 }
230 });
231 return m;
232 }
233
234
235 @Test
236 void testRank() {
237 double[][] d = { { 1, 1, 1 }, { 0, 0, 0 }, { 1, 2, 3 } };
238 RealMatrix m = new Array2DRowRealMatrix(d);
239 RRQRDecomposition qr = new RRQRDecomposition(m);
240 assertEquals(2, qr.getRank(1.0e-16));
241 }
242
243 }