Java > Open Source Codes > org > apache > commons > math > linear > BigMatrixImpl


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

A to Z: JavaDoc & Examples

---

Daily Java News & Articles

---

Open Source Projects

---

Open Source Codes

---

Free Computer Books

---

Remove Frame

---

Popular Tags