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.LocalizedCoreFormats;
26 import org.hipparchus.exception.MathIllegalArgumentException;
27 import org.hipparchus.random.RandomDataGenerator;
28 import org.junit.jupiter.api.Test;
29
30 import java.util.Random;
31
32 import static org.junit.jupiter.api.Assertions.assertEquals;
33 import static org.junit.jupiter.api.Assertions.assertTrue;
34 import static org.junit.jupiter.api.Assertions.fail;
35
36 class HessenbergTransformerTest {
37
38 private double[][] testSquare5 = {
39 { 5, 4, 3, 2, 1 },
40 { 1, 4, 0, 3, 3 },
41 { 2, 0, 3, 0, 0 },
42 { 3, 2, 1, 2, 5 },
43 { 4, 2, 1, 4, 1 }
44 };
45
46 private double[][] testSquare3 = {
47 { 2, -1, 1 },
48 { -1, 2, 1 },
49 { 1, -1, 2 }
50 };
51
52
53
54 private double[][] testRandom = {
55 { 0.680, 0.823, -0.4440, -0.2700 },
56 { -0.211, -0.605, 0.1080, 0.0268 },
57 { 0.566, -0.330, -0.0452, 0.9040 },
58 { 0.597, 0.536, 0.2580, 0.8320 }
59 };
60
61 @Test
62 void testNonSquare() {
63 try {
64 new HessenbergTransformer(MatrixUtils.createRealMatrix(new double[3][2]));
65 fail("an exception should have been thrown");
66 } catch (MathIllegalArgumentException ime) {
67 assertEquals(LocalizedCoreFormats.NON_SQUARE_MATRIX, ime.getSpecifier());
68 }
69 }
70
71 @Test
72 void testAEqualPHPt() {
73 checkAEqualPHPt(MatrixUtils.createRealMatrix(testSquare5));
74 checkAEqualPHPt(MatrixUtils.createRealMatrix(testSquare3));
75 checkAEqualPHPt(MatrixUtils.createRealMatrix(testRandom));
76 }
77
78 @Test
79 void testPOrthogonal() {
80 checkOrthogonal(new HessenbergTransformer(MatrixUtils.createRealMatrix(testSquare5)).getP());
81 checkOrthogonal(new HessenbergTransformer(MatrixUtils.createRealMatrix(testSquare3)).getP());
82 }
83
84 @Test
85 void testPTOrthogonal() {
86 checkOrthogonal(new HessenbergTransformer(MatrixUtils.createRealMatrix(testSquare5)).getPT());
87 checkOrthogonal(new HessenbergTransformer(MatrixUtils.createRealMatrix(testSquare3)).getPT());
88 }
89
90 @Test
91 void testHessenbergForm() {
92 checkHessenbergForm(new HessenbergTransformer(MatrixUtils.createRealMatrix(testSquare5)).getH());
93 checkHessenbergForm(new HessenbergTransformer(MatrixUtils.createRealMatrix(testSquare3)).getH());
94 }
95
96 @Test
97 void testRandomData() {
98 for (int run = 0; run < 100; run++) {
99 Random r = new Random(System.currentTimeMillis());
100
101
102 int size = r.nextInt(20) + 4;
103
104 double[][] data = new double[size][size];
105 for (int i = 0; i < size; i++) {
106 for (int j = 0; j < size; j++) {
107 data[i][j] = r.nextInt(100);
108 }
109 }
110
111 RealMatrix m = MatrixUtils.createRealMatrix(data);
112 RealMatrix h = checkAEqualPHPt(m);
113 checkHessenbergForm(h);
114 }
115 }
116
117 @Test
118 void testRandomDataNormalDistribution() {
119 for (int run = 0; run < 100; run++) {
120 Random r = new Random(System.currentTimeMillis());
121 RandomDataGenerator gen = new RandomDataGenerator(100);
122
123
124 int size = r.nextInt(20) + 4;
125
126 double[][] data = new double[size][size];
127 for (int i = 0; i < size; i++) {
128 for (int j = 0; j < size; j++) {
129 data[i][j] = gen.nextNormal(0.0, r.nextDouble() * 5);
130 }
131 }
132
133 RealMatrix m = MatrixUtils.createRealMatrix(data);
134 RealMatrix h = checkAEqualPHPt(m);
135 checkHessenbergForm(h);
136 }
137 }
138
139 @Test
140 void testMatricesValues5() {
141 checkMatricesValues(testSquare5,
142 new double[][] {
143 { 1.0, 0.0, 0.0, 0.0, 0.0 },
144 { 0.0, -0.182574185835055, 0.784218758628863, 0.395029040913988, -0.442289115981669 },
145 { 0.0, -0.365148371670111, -0.337950625265477, -0.374110794088820, -0.782621974707823 },
146 { 0.0, -0.547722557505166, 0.402941130124223, -0.626468266309003, 0.381019628053472 },
147 { 0.0, -0.730296743340221, -0.329285224617644, 0.558149336547665, 0.216118545309225 }
148 },
149 new double[][] {
150 { 5.0, -3.65148371670111, 2.59962019434982, -0.237003414680848, -3.13886458663398 },
151 { -5.47722557505166, 6.9, -2.29164066120599, 0.207283564429169, 0.703858369151728 },
152 { 0.0, -4.21386600008432, 2.30555659846067, 2.74935928725112, 0.857569835914113 },
153 { 0.0, 0.0, 2.86406180891882, -1.11582249161595, 0.817995267184158 },
154 { 0.0, 0.0, 0.0, 0.683518597386085, 1.91026589315528 }
155 });
156 }
157
158 @Test
159 void testMatricesValues3() {
160 checkMatricesValues(testSquare3,
161 new double[][] {
162 { 1.0, 0.0, 0.0 },
163 { 0.0, -0.707106781186547, 0.707106781186547 },
164 { 0.0, 0.707106781186547, 0.707106781186548 },
165 },
166 new double[][] {
167 { 2.0, 1.41421356237309, 0.0 },
168 { 1.41421356237310, 2.0, -1.0 },
169 { 0.0, 1.0, 2.0 },
170 });
171 }
172
173
174
175
176
177 private RealMatrix checkAEqualPHPt(RealMatrix matrix) {
178 HessenbergTransformer transformer = new HessenbergTransformer(matrix);
179 RealMatrix p = transformer.getP();
180 RealMatrix pT = transformer.getPT();
181 RealMatrix h = transformer.getH();
182
183 RealMatrix result = p.multiply(h).multiply(pT);
184 double norm = result.subtract(matrix).getNorm1();
185 assertEquals(0, norm, 1.0e-10);
186
187 for (int i = 0; i < matrix.getRowDimension(); ++i) {
188 for (int j = 0; j < matrix.getColumnDimension(); ++j) {
189 if (i > j + 1) {
190 assertEquals(matrix.getEntry(i, j), result.getEntry(i, j), 1.0e-12);
191 }
192 }
193 }
194
195 return transformer.getH();
196 }
197
198 private void checkOrthogonal(RealMatrix m) {
199 RealMatrix mTm = m.transposeMultiply(m);
200 RealMatrix id = MatrixUtils.createRealIdentityMatrix(mTm.getRowDimension());
201 assertEquals(0, mTm.subtract(id).getNorm1(), 1.0e-14);
202 }
203
204 private void checkHessenbergForm(RealMatrix m) {
205 final int rows = m.getRowDimension();
206 final int cols = m.getColumnDimension();
207 for (int i = 0; i < rows; ++i) {
208 for (int j = 0; j < cols; ++j) {
209 if (i > j + 1) {
210 assertEquals(0, m.getEntry(i, j), 1.0e-16);
211 }
212 }
213 }
214 }
215
216 private void checkMatricesValues(double[][] matrix, double[][] pRef, double[][] hRef) {
217 HessenbergTransformer transformer =
218 new HessenbergTransformer(MatrixUtils.createRealMatrix(matrix));
219
220
221 RealMatrix p = transformer.getP();
222 assertEquals(0, p.subtract(MatrixUtils.createRealMatrix(pRef)).getNorm1(), 1.0e-14);
223
224 RealMatrix h = transformer.getH();
225 assertEquals(0, h.subtract(MatrixUtils.createRealMatrix(hRef)).getNorm1(), 1.3e-14);
226
227
228 assertTrue(p == transformer.getP());
229 assertTrue(h == transformer.getH());
230 }
231 }