001 /*
002 * Licensed to the Apache Software Foundation (ASF) under one or more
003 * contributor license agreements. See the NOTICE file distributed with
004 * this work for additional information regarding copyright ownership.
005 * The ASF licenses this file to You under the Apache License, Version 2.0
006 * (the "License"); you may not use this file except in compliance with
007 * the License. You may obtain a copy of the License at
008 *
009 * http://www.apache.org/licenses/LICENSE-2.0
010 *
011 * Unless required by applicable law or agreed to in writing, software
012 * distributed under the License is distributed on an "AS IS" BASIS,
013 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
014 * See the License for the specific language governing permissions and
015 * limitations under the License.
016 */
017
018 package org.apache.commons.math.linear;
019
020 import org.apache.commons.math.linear.InvalidMatrixException;
021 import org.apache.commons.math.linear.LUDecomposition;
022 import org.apache.commons.math.linear.LUDecompositionImpl;
023 import org.apache.commons.math.linear.MatrixUtils;
024 import org.apache.commons.math.linear.RealMatrix;
025
026 import junit.framework.Test;
027 import junit.framework.TestCase;
028 import junit.framework.TestSuite;
029
030 public class LUDecompositionImplTest extends TestCase {
031 private double[][] testData = {
032 { 1.0, 2.0, 3.0},
033 { 2.0, 5.0, 3.0},
034 { 1.0, 0.0, 8.0}
035 };
036 private double[][] testDataMinus = {
037 { -1.0, -2.0, -3.0},
038 { -2.0, -5.0, -3.0},
039 { -1.0, 0.0, -8.0}
040 };
041 private double[][] luData = {
042 { 2.0, 3.0, 3.0 },
043 { 0.0, 5.0, 7.0 },
044 { 6.0, 9.0, 8.0 }
045 };
046
047 // singular matrices
048 private double[][] singular = {
049 { 2.0, 3.0 },
050 { 2.0, 3.0 }
051 };
052 private double[][] bigSingular = {
053 { 1.0, 2.0, 3.0, 4.0 },
054 { 2.0, 5.0, 3.0, 4.0 },
055 { 7.0, 3.0, 256.0, 1930.0 },
056 { 3.0, 7.0, 6.0, 8.0 }
057 }; // 4th row = 1st + 2nd
058
059 private static final double entryTolerance = 10e-16;
060
061 private static final double normTolerance = 10e-14;
062
063 public LUDecompositionImplTest(String name) {
064 super(name);
065 }
066
067 public static Test suite() {
068 TestSuite suite = new TestSuite(LUDecompositionImplTest.class);
069 suite.setName("LUDecompositionImpl Tests");
070 return suite;
071 }
072
073 /** test dimensions */
074 public void testDimensions() {
075 RealMatrix matrix = MatrixUtils.createRealMatrix(testData);
076 LUDecomposition LU = new LUDecompositionImpl(matrix);
077 assertEquals(testData.length, LU.getL().getRowDimension());
078 assertEquals(testData.length, LU.getL().getColumnDimension());
079 assertEquals(testData.length, LU.getU().getRowDimension());
080 assertEquals(testData.length, LU.getU().getColumnDimension());
081 assertEquals(testData.length, LU.getP().getRowDimension());
082 assertEquals(testData.length, LU.getP().getColumnDimension());
083
084 }
085
086 /** test non-square matrix */
087 public void testNonSquare() {
088 try {
089 new LUDecompositionImpl(MatrixUtils.createRealMatrix(new double[3][2]));
090 } catch (InvalidMatrixException ime) {
091 // expected behavior
092 } catch (Exception e) {
093 fail("wrong exception caught");
094 }
095 }
096
097 /** test PA = LU */
098 public void testPAEqualLU() {
099 RealMatrix matrix = MatrixUtils.createRealMatrix(testData);
100 LUDecomposition lu = new LUDecompositionImpl(matrix);
101 RealMatrix l = lu.getL();
102 RealMatrix u = lu.getU();
103 RealMatrix p = lu.getP();
104 double norm = l.multiply(u).subtract(p.multiply(matrix)).getNorm();
105 assertEquals(0, norm, normTolerance);
106
107 matrix = MatrixUtils.createRealMatrix(testDataMinus);
108 lu = new LUDecompositionImpl(matrix);
109 l = lu.getL();
110 u = lu.getU();
111 p = lu.getP();
112 norm = l.multiply(u).subtract(p.multiply(matrix)).getNorm();
113 assertEquals(0, norm, normTolerance);
114
115 matrix = MatrixUtils.createRealIdentityMatrix(17);
116 lu = new LUDecompositionImpl(matrix);
117 l = lu.getL();
118 u = lu.getU();
119 p = lu.getP();
120 norm = l.multiply(u).subtract(p.multiply(matrix)).getNorm();
121 assertEquals(0, norm, normTolerance);
122
123 matrix = MatrixUtils.createRealMatrix(singular);
124 lu = new LUDecompositionImpl(matrix);
125 assertFalse(lu.getSolver().isNonSingular());
126 assertNull(lu.getL());
127 assertNull(lu.getU());
128 assertNull(lu.getP());
129
130 matrix = MatrixUtils.createRealMatrix(bigSingular);
131 lu = new LUDecompositionImpl(matrix);
132 assertFalse(lu.getSolver().isNonSingular());
133 assertNull(lu.getL());
134 assertNull(lu.getU());
135 assertNull(lu.getP());
136
137 }
138
139 /** test that L is lower triangular with unit diagonal */
140 public void testLLowerTriangular() {
141 RealMatrix matrix = MatrixUtils.createRealMatrix(testData);
142 RealMatrix l = new LUDecompositionImpl(matrix).getL();
143 for (int i = 0; i < l.getRowDimension(); i++) {
144 assertEquals(l.getEntry(i, i), 1, entryTolerance);
145 for (int j = i + 1; j < l.getColumnDimension(); j++) {
146 assertEquals(l.getEntry(i, j), 0, entryTolerance);
147 }
148 }
149 }
150
151 /** test that U is upper triangular */
152 public void testUUpperTriangular() {
153 RealMatrix matrix = MatrixUtils.createRealMatrix(testData);
154 RealMatrix u = new LUDecompositionImpl(matrix).getU();
155 for (int i = 0; i < u.getRowDimension(); i++) {
156 for (int j = 0; j < i; j++) {
157 assertEquals(u.getEntry(i, j), 0, entryTolerance);
158 }
159 }
160 }
161
162 /** test that P is a permutation matrix */
163 public void testPPermutation() {
164 RealMatrix matrix = MatrixUtils.createRealMatrix(testData);
165 RealMatrix p = new LUDecompositionImpl(matrix).getP();
166
167 RealMatrix ppT = p.multiply(p.transpose());
168 RealMatrix id = MatrixUtils.createRealIdentityMatrix(p.getRowDimension());
169 assertEquals(0, ppT.subtract(id).getNorm(), normTolerance);
170
171 for (int i = 0; i < p.getRowDimension(); i++) {
172 int zeroCount = 0;
173 int oneCount = 0;
174 int otherCount = 0;
175 for (int j = 0; j < p.getColumnDimension(); j++) {
176 final double e = p.getEntry(i, j);
177 if (e == 0) {
178 ++zeroCount;
179 } else if (e == 1) {
180 ++oneCount;
181 } else {
182 ++otherCount;
183 }
184 }
185 assertEquals(p.getColumnDimension() - 1, zeroCount);
186 assertEquals(1, oneCount);
187 assertEquals(0, otherCount);
188 }
189
190 for (int j = 0; j < p.getColumnDimension(); j++) {
191 int zeroCount = 0;
192 int oneCount = 0;
193 int otherCount = 0;
194 for (int i = 0; i < p.getRowDimension(); i++) {
195 final double e = p.getEntry(i, j);
196 if (e == 0) {
197 ++zeroCount;
198 } else if (e == 1) {
199 ++oneCount;
200 } else {
201 ++otherCount;
202 }
203 }
204 assertEquals(p.getRowDimension() - 1, zeroCount);
205 assertEquals(1, oneCount);
206 assertEquals(0, otherCount);
207 }
208
209 }
210
211
212 /** test singular */
213 public void testSingular() {
214 LUDecomposition lu =
215 new LUDecompositionImpl(MatrixUtils.createRealMatrix(testData));
216 assertTrue(lu.getSolver().isNonSingular());
217 lu = new LUDecompositionImpl(MatrixUtils.createRealMatrix(singular));
218 assertFalse(lu.getSolver().isNonSingular());
219 lu = new LUDecompositionImpl(MatrixUtils.createRealMatrix(bigSingular));
220 assertFalse(lu.getSolver().isNonSingular());
221 }
222
223 /** test matrices values */
224 public void testMatricesValues1() {
225 LUDecomposition lu =
226 new LUDecompositionImpl(MatrixUtils.createRealMatrix(testData));
227 RealMatrix lRef = MatrixUtils.createRealMatrix(new double[][] {
228 { 1.0, 0.0, 0.0 },
229 { 0.5, 1.0, 0.0 },
230 { 0.5, 0.2, 1.0 }
231 });
232 RealMatrix uRef = MatrixUtils.createRealMatrix(new double[][] {
233 { 2.0, 5.0, 3.0 },
234 { 0.0, -2.5, 6.5 },
235 { 0.0, 0.0, 0.2 }
236 });
237 RealMatrix pRef = MatrixUtils.createRealMatrix(new double[][] {
238 { 0.0, 1.0, 0.0 },
239 { 0.0, 0.0, 1.0 },
240 { 1.0, 0.0, 0.0 }
241 });
242 int[] pivotRef = { 1, 2, 0 };
243
244 // check values against known references
245 RealMatrix l = lu.getL();
246 assertEquals(0, l.subtract(lRef).getNorm(), 1.0e-13);
247 RealMatrix u = lu.getU();
248 assertEquals(0, u.subtract(uRef).getNorm(), 1.0e-13);
249 RealMatrix p = lu.getP();
250 assertEquals(0, p.subtract(pRef).getNorm(), 1.0e-13);
251 int[] pivot = lu.getPivot();
252 for (int i = 0; i < pivotRef.length; ++i) {
253 assertEquals(pivotRef[i], pivot[i]);
254 }
255
256 // check the same cached instance is returned the second time
257 assertTrue(l == lu.getL());
258 assertTrue(u == lu.getU());
259 assertTrue(p == lu.getP());
260
261 }
262
263 /** test matrices values */
264 public void testMatricesValues2() {
265 LUDecomposition lu =
266 new LUDecompositionImpl(MatrixUtils.createRealMatrix(luData));
267 RealMatrix lRef = MatrixUtils.createRealMatrix(new double[][] {
268 { 1.0, 0.0, 0.0 },
269 { 0.0, 1.0, 0.0 },
270 { 1.0 / 3.0, 0.0, 1.0 }
271 });
272 RealMatrix uRef = MatrixUtils.createRealMatrix(new double[][] {
273 { 6.0, 9.0, 8.0 },
274 { 0.0, 5.0, 7.0 },
275 { 0.0, 0.0, 1.0 / 3.0 }
276 });
277 RealMatrix pRef = MatrixUtils.createRealMatrix(new double[][] {
278 { 0.0, 0.0, 1.0 },
279 { 0.0, 1.0, 0.0 },
280 { 1.0, 0.0, 0.0 }
281 });
282 int[] pivotRef = { 2, 1, 0 };
283
284 // check values against known references
285 RealMatrix l = lu.getL();
286 assertEquals(0, l.subtract(lRef).getNorm(), 1.0e-13);
287 RealMatrix u = lu.getU();
288 assertEquals(0, u.subtract(uRef).getNorm(), 1.0e-13);
289 RealMatrix p = lu.getP();
290 assertEquals(0, p.subtract(pRef).getNorm(), 1.0e-13);
291 int[] pivot = lu.getPivot();
292 for (int i = 0; i < pivotRef.length; ++i) {
293 assertEquals(pivotRef[i], pivot[i]);
294 }
295
296 // check the same cached instance is returned the second time
297 assertTrue(l == lu.getL());
298 assertTrue(u == lu.getU());
299 assertTrue(p == lu.getP());
300
301 }
302
303 }