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:
019: import org.jscience.mathematics.number.Float64;
020:
021: /**
022: * <p> This class represents an optimized {@link Matrix matrix} implementation
023: * for {@link Float64 64 bits floating-point} numbers.</p>
024: *
025: * <p> Instances of this class can be created from {@link Float64Vector},
026: * either as rows or columns if the matrix is transposed. For example:[code]
027: * Float64Vector<Rational> column0 = Float64Vector.valueOf(...);
028: * Float64Vector<Rational> column1 = Float64Vector.valueOf(...);
029: * Float64Matrix<Rational> M = Float64Matrix.valueOf(column0, column1).transpose();
030: * [/code]</p>
031: *
032: * @author <a href="mailto:jean-marie@dautelle.com">Jean-Marie Dautelle</a>
033: * @version 3.3, January 2, 2007
034: */
035: public final class Float64Matrix extends Matrix<Float64> {
036:
037: /**
038: * Holds the number of columns n.
039: */
040: int _n;;
041:
042: /**
043: * Indicates if this matrix is transposed (the rows are then the columns).
044: */
045: boolean _transposed;
046:
047: /**
048: * Holds this matrix rows (or columns when transposed).
049: */
050: final FastTable<Float64Vector> _rows = new FastTable<Float64Vector>();
051:
052: /**
053: * Returns a dense matrix from a 2-dimensional array of <code>double</code>
054: * values. The first dimension being the row and the second being the
055: * column.
056: *
057: * @param values the array of <code>double</code> values.
058: * @return the matrix having the specified elements.
059: * @throws DimensionException if rows have different length.
060: * @see Float64Vector
061: */
062: public static Float64Matrix valueOf(double[][] values) {
063: int m = values.length;
064: int n = values[0].length;
065: Float64Matrix M = Float64Matrix.newInstance(n, false);
066: for (int i = 0; i < m; i++) {
067: Float64Vector row = Float64Vector.valueOf(values[i]);
068: if (row.getDimension() != n)
069: throw new DimensionException();
070: M._rows.add(row);
071: }
072: return M;
073: }
074:
075: /**
076: * Returns a complex matrix holding the specified row vectors
077: * (column vectors if {@link #transpose transposed}).
078: *
079: * @param rows the row vectors.
080: * @return the matrix having the specified rows.
081: * @throws DimensionException if the rows do not have the same dimension.
082: */
083: public static Float64Matrix valueOf(Float64Vector... rows) {
084: final int n = rows[0].getDimension();
085: Float64Matrix M = Float64Matrix.newInstance(n, false);
086: for (int i = 0, m = rows.length; i < m; i++) {
087: Float64Vector rowi = rows[i];
088: if (rowi.getDimension() != n)
089: throw new DimensionException(
090: "All vectors must have the same dimension.");
091: M._rows.add(rowi);
092: }
093: return M;
094: }
095:
096: /**
097: * Returns a complex matrix holding the row vectors from the specified
098: * collection (column vectors if {@link #transpose transposed}).
099: *
100: * @param rows the list of row vectors.
101: * @return the matrix having the specified rows.
102: * @throws DimensionException if the rows do not have the same dimension.
103: */
104: public static Float64Matrix valueOf(List<Float64Vector> rows) {
105: final int n = rows.get(0).getDimension();
106: Float64Matrix M = Float64Matrix.newInstance(n, false);
107: Iterator<Float64Vector> iterator = rows.iterator();
108: for (int i = 0, m = rows.size(); i < m; i++) {
109: Float64Vector rowi = iterator.next();
110: if (rowi.getDimension() != n)
111: throw new DimensionException(
112: "All vectors must have the same dimension.");
113: M._rows.add(rowi);
114: }
115: return M;
116: }
117:
118: /**
119: * Returns a complex matrix equivalent to the specified matrix.
120: *
121: * @param that the matrix to convert.
122: * @return <code>that</code> or a complex matrix holding the same elements
123: * as the specified matrix.
124: */
125: public static Float64Matrix valueOf(Matrix<Float64> that) {
126: if (that instanceof Float64Matrix)
127: return (Float64Matrix) that;
128: int n = that.getNumberOfColumns();
129: int m = that.getNumberOfRows();
130: Float64Matrix M = Float64Matrix.newInstance(n, false);
131: for (int i = 0; i < m; i++) {
132: Float64Vector rowi = Float64Vector.valueOf(that.getRow(i));
133: M._rows.add(rowi);
134: }
135: return M;
136: }
137:
138: @Override
139: public int getNumberOfRows() {
140: return _transposed ? _n : _rows.size();
141: }
142:
143: @Override
144: public int getNumberOfColumns() {
145: return _transposed ? _rows.size() : _n;
146: }
147:
148: @Override
149: public Float64 get(int i, int j) {
150: return _transposed ? _rows.get(j).get(i) : _rows.get(i).get(j);
151: }
152:
153: @Override
154: public Float64Vector getRow(int i) {
155: if (!_transposed)
156: return _rows.get(i);
157: // Else transposed.
158: int n = _rows.size();
159: int m = _n;
160: if ((i < 0) || (i >= m))
161: throw new DimensionException();
162: Float64Vector V = Float64Vector.newInstance(n);
163: for (int j = 0; j < n; j++) {
164: V.set(j, _rows.get(j).get(i).doubleValue());
165: }
166: return V;
167: }
168:
169: @Override
170: public Float64Vector getColumn(int j) {
171: if (_transposed)
172: return _rows.get(j);
173: int m = _rows.size();
174: if ((j < 0) || (j >= _n))
175: throw new DimensionException();
176: Float64Vector V = Float64Vector.newInstance(m);
177: for (int i = 0; i < m; i++) {
178: V.set(i, _rows.get(i).get(j).doubleValue());
179: }
180: return V;
181: }
182:
183: @Override
184: public Float64Vector getDiagonal() {
185: int m = this .getNumberOfRows();
186: int n = this .getNumberOfColumns();
187: int dimension = MathLib.min(m, n);
188: Float64Vector V = Float64Vector.newInstance(dimension);
189: for (int i = 0; i < dimension; i++) {
190: V.set(i, this .get(i, i).doubleValue());
191: }
192: return V;
193: }
194:
195: @Override
196: public Float64Matrix opposite() {
197: Float64Matrix M = Float64Matrix.newInstance(_n, _transposed);
198: for (int i = 0, p = _rows.size(); i < p; i++) {
199: M._rows.add(_rows.get(i).opposite());
200: }
201: return M;
202: }
203:
204: @Override
205: public Float64Matrix plus(Matrix<Float64> that) {
206: if (this .getNumberOfRows() != that.getNumberOfRows())
207: throw new DimensionException();
208: Float64Matrix M = Float64Matrix.newInstance(_n, _transposed);
209: for (int i = 0, p = _rows.size(); i < p; i++) {
210: M._rows.add(_rows.get(i).plus(
211: _transposed ? that.getColumn(i) : that.getRow(i)));
212: }
213: return M;
214: }
215:
216: @Override
217: public Float64Matrix minus(Matrix<Float64> that) { // Returns more specialized type.
218: return this .plus(that.opposite());
219: }
220:
221: @Override
222: public Float64Matrix times(Float64 k) {
223: Float64Matrix M = Float64Matrix.newInstance(_n, _transposed);
224: for (int i = 0, p = _rows.size(); i < p; i++) {
225: M._rows.add(_rows.get(i).times(k));
226: }
227: return M;
228: }
229:
230: @Override
231: public Float64Vector times(Vector<Float64> v) {
232: if (v.getDimension() != this .getNumberOfColumns())
233: throw new DimensionException();
234: final int m = this .getNumberOfRows();
235: Float64Vector V = Float64Vector.newInstance(m);
236: for (int i = 0; i < m; i++) {
237: V.set(i, this .getRow(i).times(v).doubleValue());
238: }
239: return V;
240: }
241:
242: @Override
243: public Float64Matrix times(Matrix<Float64> that) {
244: final int n = this .getNumberOfColumns();
245: final int m = this .getNumberOfRows();
246: final int p = that.getNumberOfColumns();
247: if (that.getNumberOfRows() != n)
248: throw new DimensionException();
249: // Creates a mxp matrix in transposed form (p columns vectors of size m)
250: Float64Matrix M = Float64Matrix.newInstance(m, true); // Transposed.
251: M._rows.setSize(p);
252: Multiply multiply = Multiply.valueOf(this , that, 0, p, M._rows);
253: multiply.run();
254: Multiply.recycle(multiply);
255: return M;
256: }
257:
258: // Logic to multiply two matrices.
259: private static class Multiply implements Runnable {
260: private static final ObjectFactory<Multiply> FACTORY = new ObjectFactory<Multiply>() {
261:
262: @Override
263: protected Multiply create() {
264: return new Multiply();
265: }
266: };
267:
268: private Float64Matrix _left;
269:
270: private Matrix<Float64> _right;
271:
272: private int _rightColumnStart;
273:
274: private int _rightColumnEnd;
275:
276: private FastTable<Float64Vector> _columnsResult;
277:
278: static Multiply valueOf(Float64Matrix left,
279: Matrix<Float64> right, int rightColumnStart,
280: int rightColumnEnd,
281: FastTable<Float64Vector> columnsResult) {
282: Multiply multiply = Multiply.FACTORY.object();
283: multiply._left = left;
284: multiply._right = right;
285: multiply._rightColumnStart = rightColumnStart;
286: multiply._rightColumnEnd = rightColumnEnd;
287: multiply._columnsResult = columnsResult;
288: return multiply;
289: }
290:
291: static void recycle(Multiply multiply) {
292: multiply._left = null;
293: multiply._right = null;
294: multiply._columnsResult = null;
295: Multiply.FACTORY.recycle(multiply);
296: }
297:
298: public void run() {
299: if (_rightColumnEnd - _rightColumnStart < 32) { // Direct calculation.
300: FastTable<Float64Vector> rows = _left.getRows();
301: final int m = rows.size();
302: for (int j = _rightColumnStart; j < _rightColumnEnd; j++) {
303: Vector<Float64> thatColj = _right.getColumn(j);
304: Float64Vector column = Float64Vector.newInstance(m);
305: _columnsResult.set(j, column);
306: for (int i = 0; i < m; i++) {
307: column.set(i, rows.get(i).times(thatColj)
308: .doubleValue());
309: }
310: }
311: } else { // Concurrent/Recursive calculation.
312: int halfIndex = (_rightColumnStart + _rightColumnEnd) >> 1;
313: Multiply firstHalf = Multiply.valueOf(_left, _right,
314: _rightColumnStart, halfIndex, _columnsResult);
315: Multiply secondHalf = Multiply.valueOf(_left, _right,
316: halfIndex, _rightColumnEnd, _columnsResult);
317: ConcurrentContext.enter();
318: try {
319: ConcurrentContext.execute(firstHalf);
320: ConcurrentContext.execute(secondHalf);
321: } finally {
322: ConcurrentContext.exit();
323: }
324: Multiply.recycle(firstHalf);
325: Multiply.recycle(secondHalf);
326: }
327: }
328: }
329:
330: private FastTable<Float64Vector> getRows() {
331: if (!_transposed)
332: return _rows;
333: FastTable<Float64Vector> rows = FastTable.newInstance();
334: for (int i = 0; i < _n; i++) {
335: rows.add(this .getRow(i));
336: }
337: return rows;
338: }
339:
340: @Override
341: public Float64Matrix inverse() {
342: if (!isSquare())
343: throw new DimensionException("Matrix not square");
344: return Float64Matrix.valueOf(LUDecomposition.valueOf(this )
345: .inverse());
346: }
347:
348: @Override
349: public Float64 determinant() {
350: return LUDecomposition.valueOf(this ).determinant();
351: }
352:
353: @Override
354: public Float64Matrix transpose() {
355: Float64Matrix M = Float64Matrix.newInstance(_n, !_transposed);
356: M._rows.addAll(this ._rows);
357: return M;
358: }
359:
360: @Override
361: public Float64 cofactor(int i, int j) {
362: if (_transposed) {
363: int k = i;
364: i = j;
365: j = k; // Swaps i,j
366: }
367: int m = _rows.size();
368: Float64Matrix M = Float64Matrix.newInstance(m - 1, _transposed);
369: for (int k1 = 0; k1 < m; k1++) {
370: if (k1 == i)
371: continue;
372: Float64Vector row = _rows.get(k1);
373: Float64Vector V = Float64Vector.newInstance(_n - 1);
374: M._rows.add(V);
375: for (int k2 = 0, k = 0; k2 < _n; k2++) {
376: if (k2 == j)
377: continue;
378: V.set(k++, row.get(k2).doubleValue());
379: }
380: }
381: return M.determinant();
382: }
383:
384: @Override
385: public Float64Matrix adjoint() {
386: Float64Matrix M = Float64Matrix.newInstance(_n, _transposed);
387: int m = _rows.size();
388: for (int i = 0; i < m; i++) {
389: Float64Vector row = Float64Vector.newInstance(_n);
390: M._rows.add(row);
391: for (int j = 0; j < _n; j++) {
392: Float64 cofactor = _transposed ? cofactor(j, i)
393: : cofactor(i, j);
394: row.set(j, ((i + j) % 2 == 0) ? cofactor.doubleValue()
395: : cofactor.opposite().doubleValue());
396: }
397: }
398: return M.transpose();
399: }
400:
401: @Override
402: public Float64Matrix tensor(Matrix<Float64> that) {
403: return Float64Matrix.valueOf(DenseMatrix.valueOf(this ).tensor(
404: that));
405: }
406:
407: @Override
408: public Float64Vector vectorization() {
409: return Float64Vector.valueOf(DenseMatrix.valueOf(this )
410: .vectorization());
411: }
412:
413: @Override
414: public Float64Matrix copy() {
415: Float64Matrix M = Float64Matrix.newInstance(_n, _transposed);
416: for (Float64Vector row : _rows) {
417: M._rows.add(row.copy());
418: }
419: return M;
420: }
421:
422: ///////////////////////
423: // Factory creation. //
424: ///////////////////////
425:
426: static Float64Matrix newInstance(int n, boolean transposed) {
427: Float64Matrix M = FACTORY.object();
428: M._n = n;
429: M._transposed = transposed;
430: return M;
431: }
432:
433: private static ObjectFactory<Float64Matrix> FACTORY = new ObjectFactory<Float64Matrix>() {
434: @Override
435: protected Float64Matrix create() {
436: return new Float64Matrix();
437: }
438:
439: @Override
440: protected void cleanup(Float64Matrix matrix) {
441: matrix._rows.reset();
442: }
443: };
444:
445: private Float64Matrix() {
446: }
447:
448: private static final long serialVersionUID = 1L;
449:
450: }
|