001: /*
002: * JScience - Java(TM) Tools and Libraries for the Advancement of Sciences.
003: * Copyright (C) 2006 - JScience (http://jscience.org/)
004: * All rights reserved.
005: *
006: * Permission to use, copy, modify, and distribute this software is
007: * freely granted, provided that this notice is preserved.
008: */
009: package org.jscience.mathematics.vector;
010:
011: import java.util.Iterator;
012: import java.util.List;
013:
014: import javolution.context.ConcurrentContext;
015: import javolution.context.ObjectFactory;
016: import javolution.lang.MathLib;
017: import javolution.util.FastTable;
018: import org.jscience.mathematics.structure.Field;
019:
020: /**
021: * <p> This class represents a matrix made of {@link DenseVector dense
022: * vectors} (as rows). To create a dense matrix made of column vectors the
023: * {@link #transpose} method can be used.
024: * For example:[code]
025: * DenseVector<Rational> column0 = DenseVector.valueOf(...);
026: * DenseVector<Rational> column1 = DenseVector.valueOf(...);
027: * DenseMatrix<Rational> M = DenseMatrix.valueOf(column0, column1).transpose();
028: * [/code]</p>
029: * <p> As for any concrete {@link org.jscience.mathematics.structure.Structure
030: * structure}, this class is declared <code>final</code> (otherwise most
031: * operations would have to be overridden to return the appropriate type).
032: * Specialized dense matrix should sub-class {@link Matrix} directly.
033: * For example:[code]
034: * // Extension through composition.
035: * final class TriangularMatrix <F extends Field<F>> extends Matrix<F> {
036: * private DenseMatrix<F> _value; // Possible implementation.
037: * ...
038: * public TriangularMatrix opposite() { // Returns the right type.
039: * return TriangularMatrix.valueOf(_value.opposite());
040: * }
041: * ...
042: * }[/code]
043: * </p>
044: * @author <a href="mailto:jean-marie@dautelle.com">Jean-Marie Dautelle</a>
045: * @version 3.3, January 2, 2007
046: */
047: public final class DenseMatrix<F extends Field<F>> extends Matrix<F> {
048:
049: /**
050: * Holds the number of columns n.
051: */
052: int _n;
053:
054: /**
055: * Indicates if this matrix is transposed (the rows are then the columns).
056: */
057: boolean _transposed;
058:
059: /**
060: * Holds this matrix rows (or columns when transposed).
061: */
062: final FastTable<DenseVector<F>> _rows = new FastTable<DenseVector<F>>();
063:
064: /**
065: * Returns a dense matrix from the specified 2-dimensional array.
066: * The first dimension being the row and the second being the column.
067: *
068: * @param elements this matrix elements.
069: * @return the matrix having the specified elements.
070: * @throws DimensionException if rows have different length.
071: * @see DenseMatrix
072: */
073: public static <F extends Field<F>> DenseMatrix<F> valueOf(
074: F[][] elements) {
075: int m = elements.length;
076: int n = elements[0].length;
077: DenseMatrix<F> M = DenseMatrix.newInstance(n, false);
078: for (int i = 0; i < m; i++) {
079: DenseVector<F> row = DenseVector.valueOf(elements[i]);
080: if (row.getDimension() != n)
081: throw new DimensionException();
082: M._rows.add(row);
083: }
084: return M;
085: }
086:
087: /**
088: * Returns a dense matrix holding the specified row vectors
089: * (column vectors if {@link #transpose transposed}).
090: *
091: * @param rows the row vectors.
092: * @return the matrix having the specified rows.
093: * @throws DimensionException if the rows do not have the same dimension.
094: */
095: public static <F extends Field<F>> DenseMatrix<F> valueOf(
096: DenseVector<F>... rows) {
097: final int n = rows[0].getDimension();
098: DenseMatrix<F> M = DenseMatrix.newInstance(n, false);
099: for (int i = 0, m = rows.length; i < m; i++) {
100: DenseVector<F> rowi = rows[i];
101: if (rowi.getDimension() != n)
102: throw new DimensionException(
103: "All vectors must have the same dimension.");
104: M._rows.add(rowi);
105: }
106: return M;
107: }
108:
109: /**
110: * Returns a dense matrix holding the row vectors from the specified
111: * collection (column vectors if {@link #transpose transposed}).
112: *
113: * @param rows the list of row vectors.
114: * @return the matrix having the specified rows.
115: * @throws DimensionException if the rows do not have the same dimension.
116: */
117: public static <F extends Field<F>> DenseMatrix<F> valueOf(
118: List<DenseVector<F>> rows) {
119: final int n = rows.get(0).getDimension();
120: DenseMatrix<F> M = DenseMatrix.newInstance(n, false);
121: Iterator<DenseVector<F>> iterator = rows.iterator();
122: for (int i = 0, m = rows.size(); i < m; i++) {
123: DenseVector<F> rowi = iterator.next();
124: if (rowi.getDimension() != n)
125: throw new DimensionException(
126: "All vectors must have the same dimension.");
127: M._rows.add(rowi);
128: }
129: return M;
130: }
131:
132: /**
133: * Returns a dense matrix equivalent to the specified matrix.
134: *
135: * @param that the matrix to convert.
136: * @return <code>that</code> or a dense matrix holding the same elements
137: * as the specified matrix.
138: */
139: public static <F extends Field<F>> DenseMatrix<F> valueOf(
140: Matrix<F> that) {
141: if (that instanceof DenseMatrix)
142: return (DenseMatrix<F>) that;
143: int n = that.getNumberOfColumns();
144: int m = that.getNumberOfRows();
145: DenseMatrix<F> M = DenseMatrix.newInstance(n, false);
146: for (int i = 0; i < m; i++) {
147: DenseVector<F> rowi = DenseVector.valueOf(that.getRow(i));
148: M._rows.add(rowi);
149: }
150: return M;
151: }
152:
153: @Override
154: public int getNumberOfRows() {
155: return _transposed ? _n : _rows.size();
156: }
157:
158: @Override
159: public int getNumberOfColumns() {
160: return _transposed ? _rows.size() : _n;
161: }
162:
163: @Override
164: public F get(int i, int j) {
165: return _transposed ? _rows.get(j).get(i) : _rows.get(i).get(j);
166: }
167:
168: @Override
169: public DenseVector<F> getRow(int i) {
170: if (!_transposed)
171: return _rows.get(i);
172: // Else transposed.
173: int n = _rows.size();
174: int m = _n;
175: if ((i < 0) || (i >= m))
176: throw new DimensionException();
177: DenseVector<F> V = DenseVector.newInstance();
178: for (int j = 0; j < n; j++) {
179: V._elements.add(_rows.get(j).get(i));
180: }
181: return V;
182: }
183:
184: @Override
185: public DenseVector<F> getColumn(int j) {
186: if (_transposed)
187: return _rows.get(j);
188: int m = _rows.size();
189: if ((j < 0) || (j >= _n))
190: throw new DimensionException();
191: DenseVector<F> V = DenseVector.newInstance();
192: for (int i = 0; i < m; i++) {
193: V._elements.add(_rows.get(i).get(j));
194: }
195: return V;
196: }
197:
198: @Override
199: public DenseVector<F> getDiagonal() {
200: int m = this .getNumberOfRows();
201: int n = this .getNumberOfColumns();
202: int dimension = MathLib.min(m, n);
203: DenseVector<F> V = DenseVector.newInstance();
204: for (int i = 0; i < dimension; i++) {
205: V._elements.add(this .get(i, i));
206: }
207: return V;
208: }
209:
210: @Override
211: public DenseMatrix<F> opposite() {
212: DenseMatrix<F> M = DenseMatrix.newInstance(_n, _transposed);
213: for (int i = 0, p = _rows.size(); i < p; i++) {
214: M._rows.add(_rows.get(i).opposite());
215: }
216: return M;
217: }
218:
219: @Override
220: public DenseMatrix<F> plus(Matrix<F> that) {
221: if (this .getNumberOfRows() != that.getNumberOfRows())
222: throw new DimensionException();
223: DenseMatrix<F> M = DenseMatrix.newInstance(_n, _transposed);
224: for (int i = 0, p = _rows.size(); i < p; i++) {
225: M._rows.add(_rows.get(i).plus(
226: _transposed ? that.getColumn(i) : that.getRow(i)));
227: }
228: return M;
229: }
230:
231: @Override
232: public DenseMatrix<F> minus(Matrix<F> that) { // Returns more specialized type.
233: return this .plus(that.opposite());
234: }
235:
236: @Override
237: public DenseMatrix<F> times(F k) {
238: DenseMatrix<F> M = DenseMatrix.newInstance(_n, _transposed);
239: for (int i = 0, p = _rows.size(); i < p; i++) {
240: M._rows.add(_rows.get(i).times(k));
241: }
242: return M;
243: }
244:
245: @Override
246: public DenseVector<F> times(Vector<F> v) {
247: if (v.getDimension() != this .getNumberOfColumns())
248: throw new DimensionException();
249: final int m = this .getNumberOfRows();
250: DenseVector<F> V = DenseVector.newInstance();
251: for (int i = 0; i < m; i++) {
252: V._elements.add(this .getRow(i).times(v));
253: }
254: return V;
255: }
256:
257: @Override
258: public DenseMatrix<F> times(Matrix<F> that) {
259: final int n = this .getNumberOfColumns();
260: final int m = this .getNumberOfRows();
261: final int p = that.getNumberOfColumns();
262: if (that.getNumberOfRows() != n)
263: throw new DimensionException();
264: // Creates a mxp matrix in transposed form (p columns vectors of size m)
265: DenseMatrix<F> M = DenseMatrix.newInstance(m, true); // Transposed.
266: M._rows.setSize(p);
267: Multiply<F> multiply = Multiply.valueOf(this , that, 0, p,
268: M._rows);
269: multiply.run();
270: Multiply.recycle(multiply);
271: return M;
272: }
273:
274: // Logic to multiply two matrices.
275: private static class Multiply<F extends Field<F>> implements
276: Runnable {
277: private static final ObjectFactory<Multiply> FACTORY = new ObjectFactory<Multiply>() {
278:
279: @Override
280: protected Multiply create() {
281: return new Multiply();
282: }
283: };
284:
285: private DenseMatrix<F> _left;
286:
287: private Matrix<F> _right;
288:
289: private int _rightColumnStart;
290:
291: private int _rightColumnEnd;
292:
293: private FastTable<DenseVector<F>> _columnsResult;
294:
295: @SuppressWarnings("unchecked")
296: static <F extends Field<F>> Multiply<F> valueOf(
297: DenseMatrix<F> left, Matrix<F> right,
298: int rightColumnStart, int rightColumnEnd,
299: FastTable<DenseVector<F>> columnsResult) {
300: Multiply<F> multiply = Multiply.FACTORY.object();
301: multiply._left = left;
302: multiply._right = right;
303: multiply._rightColumnStart = rightColumnStart;
304: multiply._rightColumnEnd = rightColumnEnd;
305: multiply._columnsResult = columnsResult;
306: return multiply;
307: }
308:
309: static <F extends Field<F>> void recycle(Multiply<F> multiply) {
310: multiply._left = null;
311: multiply._right = null;
312: multiply._columnsResult = null;
313: Multiply.FACTORY.recycle(multiply);
314: }
315:
316: public void run() {
317: if (_rightColumnEnd - _rightColumnStart < 32) { // Direct calculation.
318: FastTable<DenseVector<F>> rows = _left.getRows();
319: final int m = rows.size();
320: for (int j = _rightColumnStart; j < _rightColumnEnd; j++) {
321: Vector<F> thatColj = _right.getColumn(j);
322: DenseVector<F> column = DenseVector.newInstance();
323: _columnsResult.set(j, column);
324: for (int i = 0; i < m; i++) {
325: column._elements.add(rows.get(i)
326: .times(thatColj));
327: }
328: }
329: } else { // Concurrent/Recursive calculation.
330: int halfIndex = (_rightColumnStart + _rightColumnEnd) >> 1;
331: Multiply<F> firstHalf = Multiply.valueOf(_left, _right,
332: _rightColumnStart, halfIndex, _columnsResult);
333: Multiply<F> secondHalf = Multiply.valueOf(_left,
334: _right, halfIndex, _rightColumnEnd,
335: _columnsResult);
336: ConcurrentContext.enter();
337: try {
338: ConcurrentContext.execute(firstHalf);
339: ConcurrentContext.execute(secondHalf);
340: } finally {
341: ConcurrentContext.exit();
342: }
343: Multiply.recycle(firstHalf);
344: Multiply.recycle(secondHalf);
345: }
346: }
347: }
348:
349: private FastTable<DenseVector<F>> getRows() {
350: if (!_transposed)
351: return _rows;
352: FastTable<DenseVector<F>> rows = FastTable.newInstance();
353: for (int i = 0; i < _n; i++) {
354: rows.add(this .getRow(i));
355: }
356: return rows;
357: }
358:
359: @Override
360: public DenseMatrix<F> inverse() {
361: if (!isSquare())
362: throw new DimensionException("Matrix not square");
363: return LUDecomposition.valueOf(this ).inverse();
364: }
365:
366: @Override
367: public F determinant() {
368: return LUDecomposition.valueOf(this ).determinant();
369: }
370:
371: @Override
372: public DenseMatrix<F> transpose() {
373: DenseMatrix<F> M = DenseMatrix.newInstance(_n, !_transposed);
374: M._rows.addAll(this ._rows);
375: return M;
376: }
377:
378: @Override
379: public F cofactor(int i, int j) {
380: if (_transposed) {
381: int k = i;
382: i = j;
383: j = k; // Swaps i,j
384: }
385: int m = _rows.size();
386: DenseMatrix<F> M = DenseMatrix.newInstance(m - 1, _transposed);
387: for (int k1 = 0; k1 < m; k1++) {
388: if (k1 == i)
389: continue;
390: DenseVector<F> row = _rows.get(k1);
391: DenseVector<F> V = DenseVector.newInstance();
392: M._rows.add(V);
393: for (int k2 = 0; k2 < _n; k2++) {
394: if (k2 == j)
395: continue;
396: V._elements.add(row.get(k2));
397: }
398: }
399: return M.determinant();
400: }
401:
402: @Override
403: public DenseMatrix<F> adjoint() {
404: DenseMatrix<F> M = DenseMatrix.newInstance(_n, _transposed);
405: int m = _rows.size();
406: for (int i = 0; i < m; i++) {
407: DenseVector<F> row = DenseVector.newInstance();
408: M._rows.add(row);
409: for (int j = 0; j < _n; j++) {
410: F cofactor = _transposed ? cofactor(j, i) : cofactor(i,
411: j);
412: row._elements.add(((i + j) % 2 == 0) ? cofactor
413: : cofactor.opposite());
414: }
415: }
416: return M.transpose();
417: }
418:
419: @Override
420: public Matrix<F> tensor(Matrix<F> that) {
421: final int this m = this .getNumberOfRows();
422: final int this n = this .getNumberOfColumns();
423: final int thatm = that.getNumberOfRows();
424: final int thatn = that.getNumberOfColumns();
425: int n = this n * thatn; // Number of columns,
426: int m = this m * thatm; // Number of rows.
427: DenseMatrix<F> M = DenseMatrix.newInstance(n, false);
428: for (int i = 0; i < m; i++) { // Row index.
429: final int i_rem_thatm = i % thatm;
430: final int i_div_thatm = i / thatm;
431: DenseVector<F> row = DenseVector.newInstance();
432: M._rows.add(row);
433: for (int j = 0; j < this n; j++) {
434: F a = this .get(i_div_thatm, j);
435: for (int k = 0; k < thatn; k++) {
436: row._elements
437: .add(a.times(that.get(i_rem_thatm, k)));
438: }
439: }
440: }
441: return M;
442: }
443:
444: @Override
445: public Vector<F> vectorization() {
446: DenseVector<F> V = DenseVector.newInstance();
447: for (int j = 0, n = this .getNumberOfColumns(); j < n; j++) {
448: Vector<F> column = this .getColumn(j);
449: for (int i = 0, m = column.getDimension(); i < m; i++) {
450: V._elements.add(column.get(i));
451: }
452: }
453: return V;
454: }
455:
456: @Override
457: public DenseMatrix<F> copy() {
458: DenseMatrix<F> M = DenseMatrix.newInstance(_n, _transposed);
459: for (DenseVector<F> row : _rows) {
460: M._rows.add(row.copy());
461: }
462: return M;
463: }
464:
465: ///////////////////////////////
466: // Package Private Utilities //
467: ///////////////////////////////
468:
469: void set(int i, int j, F e) {
470: if (_transposed) {
471: _rows.get(j)._elements.set(i, e);
472: } else {
473: _rows.get(i)._elements.set(j, e);
474: }
475: }
476:
477: ///////////////////////
478: // Factory creation. //
479: ///////////////////////
480:
481: @SuppressWarnings("unchecked")
482: static <F extends Field<F>> DenseMatrix<F> newInstance(int n,
483: boolean transposed) {
484: DenseMatrix<F> M = FACTORY.object();
485: M._n = n;
486: M._transposed = transposed;
487: return M;
488: }
489:
490: private static ObjectFactory<DenseMatrix> FACTORY = new ObjectFactory<DenseMatrix>() {
491: @Override
492: protected DenseMatrix create() {
493: return new DenseMatrix();
494: }
495:
496: @Override
497: protected void cleanup(DenseMatrix matrix) {
498: matrix._rows.reset();
499: }
500: };
501:
502: private DenseMatrix() {
503: }
504:
505: private static final long serialVersionUID = 1L;
506:
507: }
|