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 * http://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 package org.apache.commons.math.linear;
19 import java.io.Serializable;
20 import java.math.BigDecimal;
21
22 import org.apache.commons.math.MathRuntimeException;
23 import org.apache.commons.math.exception.util.LocalizedFormats;
24
25 /**
26 * Implementation of {@link BigMatrix} using a BigDecimal[][] array to store entries
27 * and <a href="http://www.math.gatech.edu/~bourbaki/math2601/Web-notes/2num.pdf">
28 * LU decompostion</a> to support linear system
29 * solution and inverse.
30 * <p>
31 * The LU decompostion is performed as needed, to support the following operations: <ul>
32 * <li>solve</li>
33 * <li>isSingular</li>
34 * <li>getDeterminant</li>
35 * <li>inverse</li> </ul></p>
36 * <p>
37 * <strong>Usage notes</strong>:<br>
38 * <ul><li>
39 * The LU decomposition is stored and reused on subsequent calls. If matrix
40 * data are modified using any of the public setXxx methods, the saved
41 * decomposition is discarded. If data are modified via references to the
42 * underlying array obtained using <code>getDataRef()</code>, then the stored
43 * LU decomposition will not be discarded. In this case, you need to
44 * explicitly invoke <code>LUDecompose()</code> to recompute the decomposition
45 * before using any of the methods above.</li>
46 * <li>
47 * As specified in the {@link BigMatrix} interface, matrix element indexing
48 * is 0-based -- e.g., <code>getEntry(0, 0)</code>
49 * returns the element in the first row, first column of the matrix.</li></ul></p>
50 *
51 * @deprecated as of 2.0, replaced by {@link Array2DRowFieldMatrix} with a {@link
52 * org.apache.commons.math.util.BigReal} parameter
53 * @version $Revision: 1042376 $ $Date: 2010-12-05 16:54:55 +0100 (dim. 05 déc. 2010) $
54 */
55 @Deprecated
56 public class BigMatrixImpl implements BigMatrix, Serializable {
57
58 /** BigDecimal 0 */
59 static final BigDecimal ZERO = new BigDecimal(0);
60
61 /** BigDecimal 1 */
62 static final BigDecimal ONE = new BigDecimal(1);
63
64 /** Bound to determine effective singularity in LU decomposition */
65 private static final BigDecimal TOO_SMALL = new BigDecimal(10E-12);
66
67 /** Serialization id */
68 private static final long serialVersionUID = -1011428905656140431L;
69
70 /** Entries of the matrix */
71 protected BigDecimal data[][] = null;
72
73 /** Entries of cached LU decomposition.
74 * All updates to data (other than luDecompose()) *must* set this to null
75 */
76 protected BigDecimal lu[][] = null;
77
78 /** Permutation associated with LU decomposition */
79 protected int[] permutation = null;
80
81 /** Parity of the permutation associated with the LU decomposition */
82 protected int parity = 1;
83
84 /** Rounding mode for divisions **/
85 private int roundingMode = BigDecimal.ROUND_HALF_UP;
86
87 /*** BigDecimal scale ***/
88 private int scale = 64;
89
90 /**
91 * Creates a matrix with no data
92 */
93 public BigMatrixImpl() {
94 }
95
96 /**
97 * Create a new BigMatrix with the supplied row and column dimensions.
98 *
99 * @param rowDimension the number of rows in the new matrix
100 * @param columnDimension the number of columns in the new matrix
101 * @throws IllegalArgumentException if row or column dimension is not
102 * positive
103 */
104 public BigMatrixImpl(int rowDimension, int columnDimension) {
105 if (rowDimension < 1 ) {
106 throw MathRuntimeException.createIllegalArgumentException(
107 LocalizedFormats.INSUFFICIENT_DIMENSION, rowDimension, 1);
108 }
109 if (columnDimension < 1) {
110 throw MathRuntimeException.createIllegalArgumentException(
111 LocalizedFormats.INSUFFICIENT_DIMENSION, columnDimension, 1);
112 }
113 data = new BigDecimal[rowDimension][columnDimension];
114 lu = null;
115 }
116
117 /**
118 * Create a new BigMatrix using <code>d</code> as the underlying
119 * data array.
120 * <p>The input array is copied, not referenced. This constructor has
121 * the same effect as calling {@link #BigMatrixImpl(BigDecimal[][], boolean)}
122 * with the second argument set to <code>true</code>.</p>
123 *
124 * @param d data for new matrix
125 * @throws IllegalArgumentException if <code>d</code> is not rectangular
126 * (not all rows have the same length) or empty
127 * @throws NullPointerException if <code>d</code> is null
128 */
129 public BigMatrixImpl(BigDecimal[][] d) {
130 this.copyIn(d);
131 lu = null;
132 }
133
134 /**
135 * Create a new BigMatrix using the input array as the underlying
136 * data array.
137 * <p>If an array is built specially in order to be embedded in a
138 * BigMatrix and not used directly, the <code>copyArray</code> may be
139 * set to <code>false</code. This will prevent the copying and improve
140 * performance as no new array will be built and no data will be copied.</p>
141 * @param d data for new matrix
142 * @param copyArray if true, the input array will be copied, otherwise
143 * it will be referenced
144 * @throws IllegalArgumentException if <code>d</code> is not rectangular
145 * (not all rows have the same length) or empty
146 * @throws NullPointerException if <code>d</code> is null
147 * @see #BigMatrixImpl(BigDecimal[][])
148 */
149 public BigMatrixImpl(BigDecimal[][] d, boolean copyArray) {
150 if (copyArray) {
151 copyIn(d);
152 } else {
153 if (d == null) {
154 throw new NullPointerException();
155 }
156 final int nRows = d.length;
157 if (nRows == 0) {
158 throw MathRuntimeException.createIllegalArgumentException(LocalizedFormats.AT_LEAST_ONE_ROW);
159 }
160
161 final int nCols = d[0].length;
162 if (nCols == 0) {
163 throw MathRuntimeException.createIllegalArgumentException(LocalizedFormats.AT_LEAST_ONE_COLUMN);
164 }
165 for (int r = 1; r < nRows; r++) {
166 if (d[r].length != nCols) {
167 throw MathRuntimeException.createIllegalArgumentException(
168 LocalizedFormats.DIFFERENT_ROWS_LENGTHS,
169 nCols, d[r].length);
170 }
171 }
172 data = d;
173 }
174 lu = null;
175 }
176
177 /**
178 * Create a new BigMatrix using <code>d</code> as the underlying
179 * data array.
180 * <p>Since the underlying array will hold <code>BigDecimal</code>
181 * instances, it will be created.</p>
182 *
183 * @param d data for new matrix
184 * @throws IllegalArgumentException if <code>d</code> is not rectangular
185 * (not all rows have the same length) or empty
186 * @throws NullPointerException if <code>d</code> is null
187 */
188 public BigMatrixImpl(double[][] d) {
189 final int nRows = d.length;
190 if (nRows == 0) {
191 throw MathRuntimeException.createIllegalArgumentException(LocalizedFormats.AT_LEAST_ONE_ROW);
192 }
193
194 final int nCols = d[0].length;
195 if (nCols == 0) {
196 throw MathRuntimeException.createIllegalArgumentException(LocalizedFormats.AT_LEAST_ONE_COLUMN);
197 }
198 for (int row = 1; row < nRows; row++) {
199 if (d[row].length != nCols) {
200 throw MathRuntimeException.createIllegalArgumentException(
201 LocalizedFormats.DIFFERENT_ROWS_LENGTHS,
202 nCols, d[row].length);
203 }
204 }
205 this.copyIn(d);
206 lu = null;
207 }
208
209 /**
210 * Create a new BigMatrix using the values represented by the strings in
211 * <code>d</code> as the underlying data array.
212 *
213 * @param d data for new matrix
214 * @throws IllegalArgumentException if <code>d</code> is not rectangular
215 * (not all rows have the same length) or empty
216 * @throws NullPointerException if <code>d</code> is null
217 */
218 public BigMatrixImpl(String[][] d) {
219 final int nRows = d.length;
220 if (nRows == 0) {
221 throw MathRuntimeException.createIllegalArgumentException(LocalizedFormats.AT_LEAST_ONE_ROW);
222 }
223
224 final int nCols = d[0].length;
225 if (nCols == 0) {
226 throw MathRuntimeException.createIllegalArgumentException(LocalizedFormats.AT_LEAST_ONE_COLUMN);
227 }
228 for (int row = 1; row < nRows; row++) {
229 if (d[row].length != nCols) {
230 throw MathRuntimeException.createIllegalArgumentException(
231 LocalizedFormats.DIFFERENT_ROWS_LENGTHS,
232 nCols, d[row].length);
233 }
234 }
235 this.copyIn(d);
236 lu = null;
237 }
238
239 /**
240 * Create a new (column) BigMatrix using <code>v</code> as the
241 * data for the unique column of the <code>v.length x 1</code> matrix
242 * created.
243 * <p>
244 * The input array is copied, not referenced.</p>
245 *
246 * @param v column vector holding data for new matrix
247 */
248 public BigMatrixImpl(BigDecimal[] v) {
249 final int nRows = v.length;
250 data = new BigDecimal[nRows][1];
251 for (int row = 0; row < nRows; row++) {
252 data[row][0] = v[row];
253 }
254 }
255
256 /**
257 * Create a new BigMatrix which is a copy of this.
258 *
259 * @return the cloned matrix
260 */
261 public BigMatrix copy() {
262 return new BigMatrixImpl(this.copyOut(), false);
263 }
264
265 /**
266 * Compute the sum of this and <code>m</code>.
267 *
268 * @param m matrix to be added
269 * @return this + m
270 * @throws IllegalArgumentException if m is not the same size as this
271 */
272 public BigMatrix add(BigMatrix m) throws IllegalArgumentException {
273 try {
274 return add((BigMatrixImpl) m);
275 } catch (ClassCastException cce) {
276
277 // safety check
278 MatrixUtils.checkAdditionCompatible(this, m);
279
280 final int rowCount = getRowDimension();
281 final int columnCount = getColumnDimension();
282 final BigDecimal[][] outData = new BigDecimal[rowCount][columnCount];
283 for (int row = 0; row < rowCount; row++) {
284 final BigDecimal[] dataRow = data[row];
285 final BigDecimal[] outDataRow = outData[row];
286 for (int col = 0; col < columnCount; col++) {
287 outDataRow[col] = dataRow[col].add(m.getEntry(row, col));
288 }
289 }
290 return new BigMatrixImpl(outData, false);
291 }
292 }
293
294 /**
295 * Compute the sum of this and <code>m</code>.
296 *
297 * @param m matrix to be added
298 * @return this + m
299 * @throws IllegalArgumentException if m is not the same size as this
300 */
301 public BigMatrixImpl add(BigMatrixImpl m) throws IllegalArgumentException {
302
303 // safety check
304 MatrixUtils.checkAdditionCompatible(this, m);
305
306 final int rowCount = getRowDimension();
307 final int columnCount = getColumnDimension();
308 final BigDecimal[][] outData = new BigDecimal[rowCount][columnCount];
309 for (int row = 0; row < rowCount; row++) {
310 final BigDecimal[] dataRow = data[row];
311 final BigDecimal[] mRow = m.data[row];
312 final BigDecimal[] outDataRow = outData[row];
313 for (int col = 0; col < columnCount; col++) {
314 outDataRow[col] = dataRow[col].add(mRow[col]);
315 }
316 }
317 return new BigMatrixImpl(outData, false);
318 }
319
320 /**
321 * Compute this minus <code>m</code>.
322 *
323 * @param m matrix to be subtracted
324 * @return this + m
325 * @throws IllegalArgumentException if m is not the same size as this
326 */
327 public BigMatrix subtract(BigMatrix m) throws IllegalArgumentException {
328 try {
329 return subtract((BigMatrixImpl) m);
330 } catch (ClassCastException cce) {
331
332 // safety check
333 MatrixUtils.checkSubtractionCompatible(this, m);
334
335 final int rowCount = getRowDimension();
336 final int columnCount = getColumnDimension();
337 final BigDecimal[][] outData = new BigDecimal[rowCount][columnCount];
338 for (int row = 0; row < rowCount; row++) {
339 final BigDecimal[] dataRow = data[row];
340 final BigDecimal[] outDataRow = outData[row];
341 for (int col = 0; col < columnCount; col++) {
342 outDataRow[col] = dataRow[col].subtract(getEntry(row, col));
343 }
344 }
345 return new BigMatrixImpl(outData, false);
346 }
347 }
348
349 /**
350 * Compute this minus <code>m</code>.
351 *
352 * @param m matrix to be subtracted
353 * @return this + m
354 * @throws IllegalArgumentException if m is not the same size as this
355 */
356 public BigMatrixImpl subtract(BigMatrixImpl m) throws IllegalArgumentException {
357
358 // safety check
359 MatrixUtils.checkSubtractionCompatible(this, m);
360
361 final int rowCount = getRowDimension();
362 final int columnCount = getColumnDimension();
363 final BigDecimal[][] outData = new BigDecimal[rowCount][columnCount];
364 for (int row = 0; row < rowCount; row++) {
365 final BigDecimal[] dataRow = data[row];
366 final BigDecimal[] mRow = m.data[row];
367 final BigDecimal[] outDataRow = outData[row];
368 for (int col = 0; col < columnCount; col++) {
369 outDataRow[col] = dataRow[col].subtract(mRow[col]);
370 }
371 }
372 return new BigMatrixImpl(outData, false);
373 }
374
375 /**
376 * Returns the result of adding d to each entry of this.
377 *
378 * @param d value to be added to each entry
379 * @return d + this
380 */
381 public BigMatrix scalarAdd(BigDecimal d) {
382 final int rowCount = getRowDimension();
383 final int columnCount = getColumnDimension();
384 final BigDecimal[][] outData = new BigDecimal[rowCount][columnCount];
385 for (int row = 0; row < rowCount; row++) {
386 final BigDecimal[] dataRow = data[row];
387 final BigDecimal[] outDataRow = outData[row];
388 for (int col = 0; col < columnCount; col++) {
389 outDataRow[col] = dataRow[col].add(d);
390 }
391 }
392 return new BigMatrixImpl(outData, false);
393 }
394
395 /**
396 * Returns the result of multiplying each entry of this by <code>d</code>
397 * @param d value to multiply all entries by
398 * @return d * this
399 */
400 public BigMatrix scalarMultiply(BigDecimal d) {
401 final int rowCount = getRowDimension();
402 final int columnCount = getColumnDimension();
403 final BigDecimal[][] outData = new BigDecimal[rowCount][columnCount];
404 for (int row = 0; row < rowCount; row++) {
405 final BigDecimal[] dataRow = data[row];
406 final BigDecimal[] outDataRow = outData[row];
407 for (int col = 0; col < columnCount; col++) {
408 outDataRow[col] = dataRow[col].multiply(d);
409 }
410 }
411 return new BigMatrixImpl(outData, false);
412 }
413
414 /**
415 * Returns the result of postmultiplying this by <code>m</code>.
416 * @param m matrix to postmultiply by
417 * @return this*m
418 * @throws IllegalArgumentException
419 * if columnDimension(this) != rowDimension(m)
420 */
421 public BigMatrix multiply(BigMatrix m) throws IllegalArgumentException {
422 try {
423 return multiply((BigMatrixImpl) m);
424 } catch (ClassCastException cce) {
425
426 // safety check
427 MatrixUtils.checkMultiplicationCompatible(this, m);
428
429 final int nRows = this.getRowDimension();
430 final int nCols = m.getColumnDimension();
431 final int nSum = this.getColumnDimension();
432 final BigDecimal[][] outData = new BigDecimal[nRows][nCols];
433 for (int row = 0; row < nRows; row++) {
434 final BigDecimal[] dataRow = data[row];
435 final BigDecimal[] outDataRow = outData[row];
436 for (int col = 0; col < nCols; col++) {
437 BigDecimal sum = ZERO;
438 for (int i = 0; i < nSum; i++) {
439 sum = sum.add(dataRow[i].multiply(m.getEntry(i, col)));
440 }
441 outDataRow[col] = sum;
442 }
443 }
444 return new BigMatrixImpl(outData, false);
445 }
446 }
447
448 /**
449 * Returns the result of postmultiplying this by <code>m</code>.
450 * @param m matrix to postmultiply by
451 * @return this*m
452 * @throws IllegalArgumentException
453 * if columnDimension(this) != rowDimension(m)
454 */
455 public BigMatrixImpl multiply(BigMatrixImpl m) throws IllegalArgumentException {
456
457 // safety check
458 MatrixUtils.checkMultiplicationCompatible(this, m);
459
460 final int nRows = this.getRowDimension();
461 final int nCols = m.getColumnDimension();
462 final int nSum = this.getColumnDimension();
463 final BigDecimal[][] outData = new BigDecimal[nRows][nCols];
464 for (int row = 0; row < nRows; row++) {
465 final BigDecimal[] dataRow = data[row];
466 final BigDecimal[] outDataRow = outData[row];
467 for (int col = 0; col < nCols; col++) {
468 BigDecimal sum = ZERO;
469 for (int i = 0; i < nSum; i++) {
470 sum = sum.add(dataRow[i].multiply(m.data[i][col]));
471 }
472 outDataRow[col] = sum;
473 }
474 }
475 return new BigMatrixImpl(outData, false);
476 }
477
478 /**
479 * Returns the result premultiplying this by <code>m</code>.
480 * @param m matrix to premultiply by
481 * @return m * this
482 * @throws IllegalArgumentException
483 * if rowDimension(this) != columnDimension(m)
484 */
485 public BigMatrix preMultiply(BigMatrix m) throws IllegalArgumentException {
486 return m.multiply(this);
487 }
488
489 /**
490 * Returns matrix entries as a two-dimensional array.
491 * <p>
492 * Makes a fresh copy of the underlying data.</p>
493 *
494 * @return 2-dimensional array of entries
495 */
496 public BigDecimal[][] getData() {
497 return copyOut();
498 }
499
500 /**
501 * Returns matrix entries as a two-dimensional array.
502 * <p>
503 * Makes a fresh copy of the underlying data converted to
504 * <code>double</code> values.</p>
505 *
506 * @return 2-dimensional array of entries
507 */
508 public double[][] getDataAsDoubleArray() {
509 final int nRows = getRowDimension();
510 final int nCols = getColumnDimension();
511 final double d[][] = new double[nRows][nCols];
512 for (int i = 0; i < nRows; i++) {
513 for (int j = 0; j < nCols; j++) {
514 d[i][j] = data[i][j].doubleValue();
515 }
516 }
517 return d;
518 }
519
520 /**
521 * Returns a reference to the underlying data array.
522 * <p>
523 * Does not make a fresh copy of the underlying data.</p>
524 *
525 * @return 2-dimensional array of entries
526 */
527 public BigDecimal[][] getDataRef() {
528 return data;
529 }
530
531 /***
532 * Gets the rounding mode for division operations
533 * The default is {@link java.math.BigDecimal#ROUND_HALF_UP}
534 * @see BigDecimal
535 * @return the rounding mode.
536 */
537 public int getRoundingMode() {
538 return roundingMode;
539 }
540
541 /***
542 * Sets the rounding mode for decimal divisions.
543 * @see BigDecimal
544 * @param roundingMode rounding mode for decimal divisions
545 */
546 public void setRoundingMode(int roundingMode) {
547 this.roundingMode = roundingMode;
548 }
549
550 /***
551 * Sets the scale for division operations.
552 * The default is 64
553 * @see BigDecimal
554 * @return the scale
555 */
556 public int getScale() {
557 return scale;
558 }
559
560 /***
561 * Sets the scale for division operations.
562 * @see BigDecimal
563 * @param scale scale for division operations
564 */
565 public void setScale(int scale) {
566 this.scale = scale;
567 }
568
569 /**
570 * Returns the <a href="http://mathworld.wolfram.com/MaximumAbsoluteRowSumNorm.html">
571 * maximum absolute row sum norm</a> of the matrix.
572 *
573 * @return norm
574 */
575 public BigDecimal getNorm() {
576 BigDecimal maxColSum = ZERO;
577 for (int col = 0; col < this.getColumnDimension(); col++) {
578 BigDecimal sum = ZERO;
579 for (int row = 0; row < this.getRowDimension(); row++) {
580 sum = sum.add(data[row][col].abs());
581 }
582 maxColSum = maxColSum.max(sum);
583 }
584 return maxColSum;
585 }
586
587 /**
588 * Gets a submatrix. Rows and columns are indicated
589 * counting from 0 to n-1.
590 *
591 * @param startRow Initial row index
592 * @param endRow Final row index
593 * @param startColumn Initial column index
594 * @param endColumn Final column index
595 * @return The subMatrix containing the data of the
596 * specified rows and columns
597 * @exception MatrixIndexException if row or column selections are not valid
598 */
599 public BigMatrix getSubMatrix(int startRow, int endRow,
600 int startColumn, int endColumn)
601 throws MatrixIndexException {
602
603 MatrixUtils.checkRowIndex(this, startRow);
604 MatrixUtils.checkRowIndex(this, endRow);
605 if (startRow > endRow) {
606 throw new MatrixIndexException(LocalizedFormats.INITIAL_ROW_AFTER_FINAL_ROW,
607 startRow, endRow);
608 }
609
610 MatrixUtils.checkColumnIndex(this, startColumn);
611 MatrixUtils.checkColumnIndex(this, endColumn);
612 if (startColumn > endColumn) {
613 throw new MatrixIndexException(LocalizedFormats.INITIAL_COLUMN_AFTER_FINAL_COLUMN,
614 startColumn, endColumn);
615 }
616
617 final BigDecimal[][] subMatrixData =
618 new BigDecimal[endRow - startRow + 1][endColumn - startColumn + 1];
619 for (int i = startRow; i <= endRow; i++) {
620 System.arraycopy(data[i], startColumn,
621 subMatrixData[i - startRow], 0,
622 endColumn - startColumn + 1);
623 }
624
625 return new BigMatrixImpl(subMatrixData, false);
626
627 }
628
629 /**
630 * Gets a submatrix. Rows and columns are indicated
631 * counting from 0 to n-1.
632 *
633 * @param selectedRows Array of row indices must be non-empty
634 * @param selectedColumns Array of column indices must be non-empty
635 * @return The subMatrix containing the data in the
636 * specified rows and columns
637 * @exception MatrixIndexException if supplied row or column index arrays
638 * are not valid
639 */
640 public BigMatrix getSubMatrix(int[] selectedRows, int[] selectedColumns)
641 throws MatrixIndexException {
642
643 if (selectedRows.length * selectedColumns.length == 0) {
644 if (selectedRows.length == 0) {
645 throw new MatrixIndexException(LocalizedFormats.EMPTY_SELECTED_ROW_INDEX_ARRAY);
646 }
647 throw new MatrixIndexException(LocalizedFormats.EMPTY_SELECTED_COLUMN_INDEX_ARRAY);
648 }
649
650 final BigDecimal[][] subMatrixData =
651 new BigDecimal[selectedRows.length][selectedColumns.length];
652 try {
653 for (int i = 0; i < selectedRows.length; i++) {
654 final BigDecimal[] subI = subMatrixData[i];
655 final BigDecimal[] dataSelectedI = data[selectedRows[i]];
656 for (int j = 0; j < selectedColumns.length; j++) {
657 subI[j] = dataSelectedI[selectedColumns[j]];
658 }
659 }
660 } catch (ArrayIndexOutOfBoundsException e) {
661 // we redo the loop with checks enabled
662 // in order to generate an appropriate message
663 for (final int row : selectedRows) {
664 MatrixUtils.checkRowIndex(this, row);
665 }
666 for (final int column : selectedColumns) {
667 MatrixUtils.checkColumnIndex(this, column);
668 }
669 }
670 return new BigMatrixImpl(subMatrixData, false);
671 }
672
673 /**
674 * Replace the submatrix starting at <code>row, column</code> using data in
675 * the input <code>subMatrix</code> array. Indexes are 0-based.
676 * <p>
677 * Example:<br>
678 * Starting with <pre>
679 * 1 2 3 4
680 * 5 6 7 8
681 * 9 0 1 2
682 * </pre>
683 * and <code>subMatrix = {{3, 4} {5,6}}</code>, invoking
684 * <code>setSubMatrix(subMatrix,1,1))</code> will result in <pre>
685 * 1 2 3 4
686 * 5 3 4 8
687 * 9 5 6 2
688 * </pre></p>
689 *
690 * @param subMatrix array containing the submatrix replacement data
691 * @param row row coordinate of the top, left element to be replaced
692 * @param column column coordinate of the top, left element to be replaced
693 * @throws MatrixIndexException if subMatrix does not fit into this
694 * matrix from element in (row, column)
695 * @throws IllegalArgumentException if <code>subMatrix</code> is not rectangular
696 * (not all rows have the same length) or empty
697 * @throws NullPointerException if <code>subMatrix</code> is null
698 * @since 1.1
699 */
700 public void setSubMatrix(BigDecimal[][] subMatrix, int row, int column)
701 throws MatrixIndexException {
702
703 final int nRows = subMatrix.length;
704 if (nRows == 0) {
705 throw MathRuntimeException.createIllegalArgumentException(LocalizedFormats.AT_LEAST_ONE_ROW);
706 }
707
708 final int nCols = subMatrix[0].length;
709 if (nCols == 0) {
710 throw MathRuntimeException.createIllegalArgumentException(LocalizedFormats.AT_LEAST_ONE_COLUMN);
711 }
712
713 for (int r = 1; r < nRows; r++) {
714 if (subMatrix[r].length != nCols) {
715 throw MathRuntimeException.createIllegalArgumentException(
716 LocalizedFormats.DIFFERENT_ROWS_LENGTHS,
717 nCols, subMatrix[r].length);
718 }
719 }
720
721 if (data == null) {
722 if (row > 0) {
723 throw MathRuntimeException.createIllegalStateException(
724 LocalizedFormats.FIRST_ROWS_NOT_INITIALIZED_YET,
725 row);
726 }
727 if (column > 0) {
728 throw MathRuntimeException.createIllegalStateException(
729 LocalizedFormats.FIRST_COLUMNS_NOT_INITIALIZED_YET,
730 column);
731 }
732 data = new BigDecimal[nRows][nCols];
733 System.arraycopy(subMatrix, 0, data, 0, subMatrix.length);
734 } else {
735 MatrixUtils.checkRowIndex(this, row);
736 MatrixUtils.checkColumnIndex(this, column);
737 MatrixUtils.checkRowIndex(this, nRows + row - 1);
738 MatrixUtils.checkColumnIndex(this, nCols + column - 1);
739 }
740 for (int i = 0; i < nRows; i++) {
741 System.arraycopy(subMatrix[i], 0, data[row + i], column, nCols);
742 }
743
744 lu = null;
745
746 }
747
748 /**
749 * Returns the entries in row number <code>row</code>
750 * as a row matrix. Row indices start at 0.
751 *
752 * @param row the row to be fetched
753 * @return row matrix
754 * @throws MatrixIndexException if the specified row index is invalid
755 */
756 public BigMatrix getRowMatrix(int row) throws MatrixIndexException {
757 MatrixUtils.checkRowIndex(this, row);
758 final int ncols = this.getColumnDimension();
759 final BigDecimal[][] out = new BigDecimal[1][ncols];
760 System.arraycopy(data[row], 0, out[0], 0, ncols);
761 return new BigMatrixImpl(out, false);
762 }
763
764 /**
765 * Returns the entries in column number <code>column</code>
766 * as a column matrix. Column indices start at 0.
767 *
768 * @param column the column to be fetched
769 * @return column matrix
770 * @throws MatrixIndexException if the specified column index is invalid
771 */
772 public BigMatrix getColumnMatrix(int column) throws MatrixIndexException {
773 MatrixUtils.checkColumnIndex(this, column);
774 final int nRows = this.getRowDimension();
775 final BigDecimal[][] out = new BigDecimal[nRows][1];
776 for (int row = 0; row < nRows; row++) {
777 out[row][0] = data[row][column];
778 }
779 return new BigMatrixImpl(out, false);
780 }
781
782 /**
783 * Returns the entries in row number <code>row</code> as an array.
784 * <p>
785 * Row indices start at 0. A <code>MatrixIndexException</code> is thrown
786 * unless <code>0 <= row < rowDimension.</code></p>
787 *
788 * @param row the row to be fetched
789 * @return array of entries in the row
790 * @throws MatrixIndexException if the specified row index is not valid
791 */
792 public BigDecimal[] getRow(int row) throws MatrixIndexException {
793 MatrixUtils.checkRowIndex(this, row);
794 final int ncols = this.getColumnDimension();
795 final BigDecimal[] out = new BigDecimal[ncols];
796 System.arraycopy(data[row], 0, out, 0, ncols);
797 return out;
798 }
799
800 /**
801 * Returns the entries in row number <code>row</code> as an array
802 * of double values.
803 * <p>
804 * Row indices start at 0. A <code>MatrixIndexException</code> is thrown
805 * unless <code>0 <= row < rowDimension.</code></p>
806 *
807 * @param row the row to be fetched
808 * @return array of entries in the row
809 * @throws MatrixIndexException if the specified row index is not valid
810 */
811 public double[] getRowAsDoubleArray(int row) throws MatrixIndexException {
812 MatrixUtils.checkRowIndex(this, row);
813 final int ncols = this.getColumnDimension();
814 final double[] out = new double[ncols];
815 for (int i=0;i<ncols;i++) {
816 out[i] = data[row][i].doubleValue();
817 }
818 return out;
819 }
820
821 /**
822 * Returns the entries in column number <code>col</code> as an array.
823 * <p>
824 * Column indices start at 0. A <code>MatrixIndexException</code> is thrown
825 * unless <code>0 <= column < columnDimension.</code></p>
826 *
827 * @param col the column to be fetched
828 * @return array of entries in the column
829 * @throws MatrixIndexException if the specified column index is not valid
830 */
831 public BigDecimal[] getColumn(int col) throws MatrixIndexException {
832 MatrixUtils.checkColumnIndex(this, col);
833 final int nRows = this.getRowDimension();
834 final BigDecimal[] out = new BigDecimal[nRows];
835 for (int i = 0; i < nRows; i++) {
836 out[i] = data[i][col];
837 }
838 return out;
839 }
840
841 /**
842 * Returns the entries in column number <code>col</code> as an array
843 * of double values.
844 * <p>
845 * Column indices start at 0. A <code>MatrixIndexException</code> is thrown
846 * unless <code>0 <= column < columnDimension.</code></p>
847 *
848 * @param col the column to be fetched
849 * @return array of entries in the column
850 * @throws MatrixIndexException if the specified column index is not valid
851 */
852 public double[] getColumnAsDoubleArray(int col) throws MatrixIndexException {
853 MatrixUtils.checkColumnIndex(this, col);
854 final int nrows = this.getRowDimension();
855 final double[] out = new double[nrows];
856 for (int i=0;i<nrows;i++) {
857 out[i] = data[i][col].doubleValue();
858 }
859 return out;
860 }
861
862 /**
863 * Returns the entry in the specified row and column.
864 * <p>
865 * Row and column indices start at 0 and must satisfy
866 * <ul>
867 * <li><code>0 <= row < rowDimension</code></li>
868 * <li><code> 0 <= column < columnDimension</code></li>
869 * </ul>
870 * otherwise a <code>MatrixIndexException</code> is thrown.</p>
871 *
872 * @param row row location of entry to be fetched
873 * @param column column location of entry to be fetched
874 * @return matrix entry in row,column
875 * @throws MatrixIndexException if the row or column index is not valid
876 */
877 public BigDecimal getEntry(int row, int column)
878 throws MatrixIndexException {
879 try {
880 return data[row][column];
881 } catch (ArrayIndexOutOfBoundsException e) {
882 throw new MatrixIndexException(
883 LocalizedFormats.NO_SUCH_MATRIX_ENTRY,
884 row, column, getRowDimension(), getColumnDimension());
885 }
886 }
887
888 /**
889 * Returns the entry in the specified row and column as a double.
890 * <p>
891 * Row and column indices start at 0 and must satisfy
892 * <ul>
893 * <li><code>0 <= row < rowDimension</code></li>
894 * <li><code> 0 <= column < columnDimension</code></li>
895 * </ul>
896 * otherwise a <code>MatrixIndexException</code> is thrown.</p>
897 *
898 * @param row row location of entry to be fetched
899 * @param column column location of entry to be fetched
900 * @return matrix entry in row,column
901 * @throws MatrixIndexException if the row
902 * or column index is not valid
903 */
904 public double getEntryAsDouble(int row, int column) throws MatrixIndexException {
905 return getEntry(row,column).doubleValue();
906 }
907
908 /**
909 * Returns the transpose matrix.
910 *
911 * @return transpose matrix
912 */
913 public BigMatrix transpose() {
914 final int nRows = this.getRowDimension();
915 final int nCols = this.getColumnDimension();
916 final BigDecimal[][] outData = new BigDecimal[nCols][nRows];
917 for (int row = 0; row < nRows; row++) {
918 final BigDecimal[] dataRow = data[row];
919 for (int col = 0; col < nCols; col++) {
920 outData[col][row] = dataRow[col];
921 }
922 }
923 return new BigMatrixImpl(outData, false);
924 }
925
926 /**
927 * Returns the inverse matrix if this matrix is invertible.
928 *
929 * @return inverse matrix
930 * @throws InvalidMatrixException if this is not invertible
931 */
932 public BigMatrix inverse() throws InvalidMatrixException {
933 return solve(MatrixUtils.createBigIdentityMatrix(getRowDimension()));
934 }
935
936 /**
937 * Returns the determinant of this matrix.
938 *
939 * @return determinant
940 * @throws InvalidMatrixException if matrix is not square
941 */
942 public BigDecimal getDeterminant() throws InvalidMatrixException {
943 if (!isSquare()) {
944 throw new NonSquareMatrixException(getRowDimension(), getColumnDimension());
945 }
946 if (isSingular()) { // note: this has side effect of attempting LU decomp if lu == null
947 return ZERO;
948 } else {
949 BigDecimal det = (parity == 1) ? ONE : ONE.negate();
950 for (int i = 0; i < getRowDimension(); i++) {
951 det = det.multiply(lu[i][i]);
952 }
953 return det;
954 }
955 }
956
957 /**
958 * Is this a square matrix?
959 * @return true if the matrix is square (rowDimension = columnDimension)
960 */
961 public boolean isSquare() {
962 return getColumnDimension() == getRowDimension();
963 }
964
965 /**
966 * Is this a singular matrix?
967 * @return true if the matrix is singular
968 */
969 public boolean isSingular() {
970 if (lu == null) {
971 try {
972 luDecompose();
973 return false;
974 } catch (InvalidMatrixException ex) {
975 return true;
976 }
977 } else { // LU decomp must have been successfully performed
978 return false; // so the matrix is not singular
979 }
980 }
981
982 /**
983 * Returns the number of rows in the matrix.
984 *
985 * @return rowDimension
986 */
987 public int getRowDimension() {
988 return data.length;
989 }
990
991 /**
992 * Returns the number of columns in the matrix.
993 *
994 * @return columnDimension
995 */
996 public int getColumnDimension() {
997 return data[0].length;
998 }
999
1000 /**
1001 * Returns the <a href="http://mathworld.wolfram.com/MatrixTrace.html">
1002 * trace</a> of the matrix (the sum of the elements on the main diagonal).
1003 *
1004 * @return trace
1005 *
1006 * @throws IllegalArgumentException if this matrix is not square.
1007 */
1008 public BigDecimal getTrace() throws IllegalArgumentException {
1009 if (!isSquare()) {
1010 throw new NonSquareMatrixException(getRowDimension(), getColumnDimension());
1011 }
1012 BigDecimal trace = data[0][0];
1013 for (int i = 1; i < this.getRowDimension(); i++) {
1014 trace = trace.add(data[i][i]);
1015 }
1016 return trace;
1017 }
1018
1019 /**
1020 * Returns the result of multiplying this by the vector <code>v</code>.
1021 *
1022 * @param v the vector to operate on
1023 * @return this*v
1024 * @throws IllegalArgumentException if columnDimension != v.size()
1025 */
1026 public BigDecimal[] operate(BigDecimal[] v) throws IllegalArgumentException {
1027 if (v.length != getColumnDimension()) {
1028 throw MathRuntimeException.createIllegalArgumentException(
1029 LocalizedFormats.VECTOR_LENGTH_MISMATCH,
1030 v.length, getColumnDimension() );
1031 }
1032 final int nRows = this.getRowDimension();
1033 final int nCols = this.getColumnDimension();
1034 final BigDecimal[] out = new BigDecimal[nRows];
1035 for (int row = 0; row < nRows; row++) {
1036 BigDecimal sum = ZERO;
1037 for (int i = 0; i < nCols; i++) {
1038 sum = sum.add(data[row][i].multiply(v[i]));
1039 }
1040 out[row] = sum;
1041 }
1042 return out;
1043 }
1044
1045 /**
1046 * Returns the result of multiplying this by the vector <code>v</code>.
1047 *
1048 * @param v the vector to operate on
1049 * @return this*v
1050 * @throws IllegalArgumentException if columnDimension != v.size()
1051 */
1052 public BigDecimal[] operate(double[] v) throws IllegalArgumentException {
1053 final BigDecimal bd[] = new BigDecimal[v.length];
1054 for (int i = 0; i < bd.length; i++) {
1055 bd[i] = new BigDecimal(v[i]);
1056 }
1057 return operate(bd);
1058 }
1059
1060 /**
1061 * Returns the (row) vector result of premultiplying this by the vector <code>v</code>.
1062 *
1063 * @param v the row vector to premultiply by
1064 * @return v*this
1065 * @throws IllegalArgumentException if rowDimension != v.size()
1066 */
1067 public BigDecimal[] preMultiply(BigDecimal[] v) throws IllegalArgumentException {
1068 final int nRows = this.getRowDimension();
1069 if (v.length != nRows) {
1070 throw MathRuntimeException.createIllegalArgumentException(
1071 LocalizedFormats.VECTOR_LENGTH_MISMATCH,
1072 v.length, nRows );
1073 }
1074 final int nCols = this.getColumnDimension();
1075 final BigDecimal[] out = new BigDecimal[nCols];
1076 for (int col = 0; col < nCols; col++) {
1077 BigDecimal sum = ZERO;
1078 for (int i = 0; i < nRows; i++) {
1079 sum = sum.add(data[i][col].multiply(v[i]));
1080 }
1081 out[col] = sum;
1082 }
1083 return out;
1084 }
1085
1086 /**
1087 * Returns a matrix of (column) solution vectors for linear systems with
1088 * coefficient matrix = this and constant vectors = columns of
1089 * <code>b</code>.
1090 *
1091 * @param b array of constants forming RHS of linear systems to
1092 * to solve
1093 * @return solution array
1094 * @throws IllegalArgumentException if this.rowDimension != row dimension
1095 * @throws InvalidMatrixException if this matrix is not square or is singular
1096 */
1097 public BigDecimal[] solve(BigDecimal[] b) throws IllegalArgumentException, InvalidMatrixException {
1098 final int nRows = this.getRowDimension();
1099 if (b.length != nRows) {
1100 throw MathRuntimeException.createIllegalArgumentException(
1101 LocalizedFormats.VECTOR_LENGTH_MISMATCH,
1102 b.length, nRows);
1103 }
1104 final BigMatrix bMatrix = new BigMatrixImpl(b);
1105 final BigDecimal[][] solution = ((BigMatrixImpl) (solve(bMatrix))).getDataRef();
1106 final BigDecimal[] out = new BigDecimal[nRows];
1107 for (int row = 0; row < nRows; row++) {
1108 out[row] = solution[row][0];
1109 }
1110 return out;
1111 }
1112
1113 /**
1114 * Returns a matrix of (column) solution vectors for linear systems with
1115 * coefficient matrix = this and constant vectors = columns of
1116 * <code>b</code>.
1117 *
1118 * @param b array of constants forming RHS of linear systems to
1119 * to solve
1120 * @return solution array
1121 * @throws IllegalArgumentException if this.rowDimension != row dimension
1122 * @throws InvalidMatrixException if this matrix is not square or is singular
1123 */
1124 public BigDecimal[] solve(double[] b) throws IllegalArgumentException, InvalidMatrixException {
1125 final BigDecimal bd[] = new BigDecimal[b.length];
1126 for (int i = 0; i < bd.length; i++) {
1127 bd[i] = new BigDecimal(b[i]);
1128 }
1129 return solve(bd);
1130 }
1131
1132 /**
1133 * Returns a matrix of (column) solution vectors for linear systems with
1134 * coefficient matrix = this and constant vectors = columns of
1135 * <code>b</code>.
1136 *
1137 * @param b matrix of constant vectors forming RHS of linear systems to
1138 * to solve
1139 * @return matrix of solution vectors
1140 * @throws IllegalArgumentException if this.rowDimension != row dimension
1141 * @throws InvalidMatrixException if this matrix is not square or is singular
1142 */
1143 public BigMatrix solve(BigMatrix b) throws IllegalArgumentException, InvalidMatrixException {
1144 if (b.getRowDimension() != getRowDimension()) {
1145 throw MathRuntimeException.createIllegalArgumentException(
1146 LocalizedFormats.DIMENSIONS_MISMATCH_2x2,
1147 b.getRowDimension(), b.getColumnDimension(), getRowDimension(), "n");
1148 }
1149 if (!isSquare()) {
1150 throw new NonSquareMatrixException(getRowDimension(), getColumnDimension());
1151 }
1152 if (this.isSingular()) { // side effect: compute LU decomp
1153 throw new SingularMatrixException();
1154 }
1155
1156 final int nCol = this.getColumnDimension();
1157 final int nColB = b.getColumnDimension();
1158 final int nRowB = b.getRowDimension();
1159
1160 // Apply permutations to b
1161 final BigDecimal[][] bp = new BigDecimal[nRowB][nColB];
1162 for (int row = 0; row < nRowB; row++) {
1163 final BigDecimal[] bpRow = bp[row];
1164 for (int col = 0; col < nColB; col++) {
1165 bpRow[col] = b.getEntry(permutation[row], col);
1166 }
1167 }
1168
1169 // Solve LY = b
1170 for (int col = 0; col < nCol; col++) {
1171 for (int i = col + 1; i < nCol; i++) {
1172 final BigDecimal[] bpI = bp[i];
1173 final BigDecimal[] luI = lu[i];
1174 for (int j = 0; j < nColB; j++) {
1175 bpI[j] = bpI[j].subtract(bp[col][j].multiply(luI[col]));
1176 }
1177 }
1178 }
1179
1180 // Solve UX = Y
1181 for (int col = nCol - 1; col >= 0; col--) {
1182 final BigDecimal[] bpCol = bp[col];
1183 final BigDecimal luDiag = lu[col][col];
1184 for (int j = 0; j < nColB; j++) {
1185 bpCol[j] = bpCol[j].divide(luDiag, scale, roundingMode);
1186 }
1187 for (int i = 0; i < col; i++) {
1188 final BigDecimal[] bpI = bp[i];
1189 final BigDecimal[] luI = lu[i];
1190 for (int j = 0; j < nColB; j++) {
1191 bpI[j] = bpI[j].subtract(bp[col][j].multiply(luI[col]));
1192 }
1193 }
1194 }
1195
1196 return new BigMatrixImpl(bp, false);
1197
1198 }
1199
1200 /**
1201 * Computes a new
1202 * <a href="http://www.math.gatech.edu/~bourbaki/math2601/Web-notes/2num.pdf">
1203 * LU decompostion</a> for this matrix, storing the result for use by other methods.
1204 * <p>
1205 * <strong>Implementation Note</strong>:<br>
1206 * Uses <a href="http://www.damtp.cam.ac.uk/user/fdl/people/sd/lectures/nummeth98/linear.htm">
1207 * Crout's algortithm</a>, with partial pivoting.</p>
1208 * <p>
1209 * <strong>Usage Note</strong>:<br>
1210 * This method should rarely be invoked directly. Its only use is
1211 * to force recomputation of the LU decomposition when changes have been
1212 * made to the underlying data using direct array references. Changes
1213 * made using setXxx methods will trigger recomputation when needed
1214 * automatically.</p>
1215 *
1216 * @throws InvalidMatrixException if the matrix is non-square or singular.
1217 */
1218 public void luDecompose() throws InvalidMatrixException {
1219
1220 final int nRows = this.getRowDimension();
1221 final int nCols = this.getColumnDimension();
1222 if (nRows != nCols) {
1223 throw new NonSquareMatrixException(getRowDimension(), getColumnDimension());
1224 }
1225 lu = this.getData();
1226
1227 // Initialize permutation array and parity
1228 permutation = new int[nRows];
1229 for (int row = 0; row < nRows; row++) {
1230 permutation[row] = row;
1231 }
1232 parity = 1;
1233
1234 // Loop over columns
1235 for (int col = 0; col < nCols; col++) {
1236
1237 BigDecimal sum = ZERO;
1238
1239 // upper
1240 for (int row = 0; row < col; row++) {
1241 final BigDecimal[] luRow = lu[row];
1242 sum = luRow[col];
1243 for (int i = 0; i < row; i++) {
1244 sum = sum.subtract(luRow[i].multiply(lu[i][col]));
1245 }
1246 luRow[col] = sum;
1247 }
1248
1249 // lower
1250 int max = col; // permutation row
1251 BigDecimal largest = ZERO;
1252 for (int row = col; row < nRows; row++) {
1253 final BigDecimal[] luRow = lu[row];
1254 sum = luRow[col];
1255 for (int i = 0; i < col; i++) {
1256 sum = sum.subtract(luRow[i].multiply(lu[i][col]));
1257 }
1258 luRow[col] = sum;
1259
1260 // maintain best permutation choice
1261 if (sum.abs().compareTo(largest) == 1) {
1262 largest = sum.abs();
1263 max = row;
1264 }
1265 }
1266
1267 // Singularity check
1268 if (lu[max][col].abs().compareTo(TOO_SMALL) <= 0) {
1269 lu = null;
1270 throw new SingularMatrixException();
1271 }
1272
1273 // Pivot if necessary
1274 if (max != col) {
1275 BigDecimal tmp = ZERO;
1276 for (int i = 0; i < nCols; i++) {
1277 tmp = lu[max][i];
1278 lu[max][i] = lu[col][i];
1279 lu[col][i] = tmp;
1280 }
1281 int temp = permutation[max];
1282 permutation[max] = permutation[col];
1283 permutation[col] = temp;
1284 parity = -parity;
1285 }
1286
1287 // Divide the lower elements by the "winning" diagonal elt.
1288 final BigDecimal luDiag = lu[col][col];
1289 for (int row = col + 1; row < nRows; row++) {
1290 final BigDecimal[] luRow = lu[row];
1291 luRow[col] = luRow[col].divide(luDiag, scale, roundingMode);
1292 }
1293
1294 }
1295
1296 }
1297
1298 /**
1299 * Get a string representation for this matrix.
1300 * @return a string representation for this matrix
1301 */
1302 @Override
1303 public String toString() {
1304 StringBuilder res = new StringBuilder();
1305 res.append("BigMatrixImpl{");
1306 if (data != null) {
1307 for (int i = 0; i < data.length; i++) {
1308 if (i > 0) {
1309 res.append(",");
1310 }
1311 res.append("{");
1312 for (int j = 0; j < data[0].length; j++) {
1313 if (j > 0) {
1314 res.append(",");
1315 }
1316 res.append(data[i][j]);
1317 }
1318 res.append("}");
1319 }
1320 }
1321 res.append("}");
1322 return res.toString();
1323 }
1324
1325 /**
1326 * Returns true iff <code>object</code> is a
1327 * <code>BigMatrixImpl</code> instance with the same dimensions as this
1328 * and all corresponding matrix entries are equal. BigDecimal.equals
1329 * is used to compare corresponding entries.
1330 *
1331 * @param object the object to test equality against.
1332 * @return true if object equals this
1333 */
1334 @Override
1335 public boolean equals(Object object) {
1336 if (object == this ) {
1337 return true;
1338 }
1339 if (object instanceof BigMatrixImpl == false) {
1340 return false;
1341 }
1342 final BigMatrix m = (BigMatrix) object;
1343 final int nRows = getRowDimension();
1344 final int nCols = getColumnDimension();
1345 if (m.getColumnDimension() != nCols || m.getRowDimension() != nRows) {
1346 return false;
1347 }
1348 for (int row = 0; row < nRows; row++) {
1349 final BigDecimal[] dataRow = data[row];
1350 for (int col = 0; col < nCols; col++) {
1351 if (!dataRow[col].equals(m.getEntry(row, col))) {
1352 return false;
1353 }
1354 }
1355 }
1356 return true;
1357 }
1358
1359 /**
1360 * Computes a hashcode for the matrix.
1361 *
1362 * @return hashcode for matrix
1363 */
1364 @Override
1365 public int hashCode() {
1366 int ret = 7;
1367 final int nRows = getRowDimension();
1368 final int nCols = getColumnDimension();
1369 ret = ret * 31 + nRows;
1370 ret = ret * 31 + nCols;
1371 for (int row = 0; row < nRows; row++) {
1372 final BigDecimal[] dataRow = data[row];
1373 for (int col = 0; col < nCols; col++) {
1374 ret = ret * 31 + (11 * (row+1) + 17 * (col+1)) *
1375 dataRow[col].hashCode();
1376 }
1377 }
1378 return ret;
1379 }
1380
1381 //------------------------ Protected methods
1382
1383 /**
1384 * Returns the LU decomposition as a BigMatrix.
1385 * Returns a fresh copy of the cached LU matrix if this has been computed;
1386 * otherwise the composition is computed and cached for use by other methods.
1387 * Since a copy is returned in either case, changes to the returned matrix do not
1388 * affect the LU decomposition property.
1389 * <p>
1390 * The matrix returned is a compact representation of the LU decomposition.
1391 * Elements below the main diagonal correspond to entries of the "L" matrix;
1392 * elements on and above the main diagonal correspond to entries of the "U"
1393 * matrix.</p>
1394 * <p>
1395 * Example: <pre>
1396 *
1397 * Returned matrix L U
1398 * 2 3 1 1 0 0 2 3 1
1399 * 5 4 6 5 1 0 0 4 6
1400 * 1 7 8 1 7 1 0 0 8
1401 * </pre>
1402 *
1403 * The L and U matrices satisfy the matrix equation LU = permuteRows(this), <br>
1404 * where permuteRows reorders the rows of the matrix to follow the order determined
1405 * by the <a href=#getPermutation()>permutation</a> property.</p>
1406 *
1407 * @return LU decomposition matrix
1408 * @throws InvalidMatrixException if the matrix is non-square or singular.
1409 */
1410 protected BigMatrix getLUMatrix() throws InvalidMatrixException {
1411 if (lu == null) {
1412 luDecompose();
1413 }
1414 return new BigMatrixImpl(lu);
1415 }
1416
1417 /**
1418 * Returns the permutation associated with the lu decomposition.
1419 * The entries of the array represent a permutation of the numbers 0, ... , nRows - 1.
1420 * <p>
1421 * Example:
1422 * permutation = [1, 2, 0] means current 2nd row is first, current third row is second
1423 * and current first row is last.</p>
1424 * <p>
1425 * Returns a fresh copy of the array.</p>
1426 *
1427 * @return the permutation
1428 */
1429 protected int[] getPermutation() {
1430 final int[] out = new int[permutation.length];
1431 System.arraycopy(permutation, 0, out, 0, permutation.length);
1432 return out;
1433 }
1434
1435 //------------------------ Private methods
1436
1437 /**
1438 * Returns a fresh copy of the underlying data array.
1439 *
1440 * @return a copy of the underlying data array.
1441 */
1442 private BigDecimal[][] copyOut() {
1443 final int nRows = this.getRowDimension();
1444 final BigDecimal[][] out = new BigDecimal[nRows][this.getColumnDimension()];
1445 // can't copy 2-d array in one shot, otherwise get row references
1446 for (int i = 0; i < nRows; i++) {
1447 System.arraycopy(data[i], 0, out[i], 0, data[i].length);
1448 }
1449 return out;
1450 }
1451
1452 /**
1453 * Replaces data with a fresh copy of the input array.
1454 * <p>
1455 * Verifies that the input array is rectangular and non-empty.</p>
1456 *
1457 * @param in data to copy in
1458 * @throws IllegalArgumentException if input array is emtpy or not
1459 * rectangular
1460 * @throws NullPointerException if input array is null
1461 */
1462 private void copyIn(BigDecimal[][] in) {
1463 setSubMatrix(in,0,0);
1464 }
1465
1466 /**
1467 * Replaces data with a fresh copy of the input array.
1468 *
1469 * @param in data to copy in
1470 */
1471 private void copyIn(double[][] in) {
1472 final int nRows = in.length;
1473 final int nCols = in[0].length;
1474 data = new BigDecimal[nRows][nCols];
1475 for (int i = 0; i < nRows; i++) {
1476 final BigDecimal[] dataI = data[i];
1477 final double[] inI = in[i];
1478 for (int j = 0; j < nCols; j++) {
1479 dataI[j] = new BigDecimal(inI[j]);
1480 }
1481 }
1482 lu = null;
1483 }
1484
1485 /**
1486 * Replaces data with BigDecimals represented by the strings in the input
1487 * array.
1488 *
1489 * @param in data to copy in
1490 */
1491 private void copyIn(String[][] in) {
1492 final int nRows = in.length;
1493 final int nCols = in[0].length;
1494 data = new BigDecimal[nRows][nCols];
1495 for (int i = 0; i < nRows; i++) {
1496 final BigDecimal[] dataI = data[i];
1497 final String[] inI = in[i];
1498 for (int j = 0; j < nCols; j++) {
1499 dataI[j] = new BigDecimal(inI[j]);
1500 }
1501 }
1502 lu = null;
1503 }
1504
1505 }